時系列解析 - ryamadaの遺伝学・遺伝統計学メモ

メモ

blog <- read.table("blog.txt",col.names=c("統計遺伝学","応用数学"))

ryamada.two <- ts(blog, start=c(2005,1), end=c(2016, 11), frequency=12)

plot(ryamada.two,plot.type="single",col=1:2)


ryamada22 <- ts(blog$統計遺伝学, start=c(2005,1), end=c(2016, 11), frequency=12)
ryamada <- ts(blog$応用数学, start=c(2005,1), end=c(2016, 11), frequency=12)

ryamadaall <- ts(blog$統計遺伝学+blog$応用数学, start=c(2005,1), end=c(2016, 12), frequency=12)

# Seasonal decomposition
fit <- stl(ryamada, s.window="period")
fit22 <- stl(ryamada22, s.window="period")
fitall <- stl(ryamadaall,s.window="period")

ts.data<- ryamada.two
ts.seasonal <- ts(cbind(fit22[[1]][,1],fit[[1]][,1]), start=c(2005,1), end=c(2016, 11), frequency=12)

ts.trend <- ts(cbind(fit22[[1]][,2],fit[[1]][,2]), start=c(2005,1), end=c(2016, 11), frequency=12)

ts.residual <- ts(cbind(fit22[[1]][,3],fit[[1]][,3]), start=c(2005,1), end=c(2016, 11), frequency=12)
op <- par(mfcol = c(2,2))
op <- par(mar = c(0, 4, 0, 3), oma = c(5, 0, 4, 0), mfcol = c(4,1))
plot(ts.data,plot.type="single",col=1:2)
plot(ts.seasonal,plot.type="single",col=1:2)
plot(ts.trend,plot.type="single",col=1:2)
plot(ts.residual,plot.type="single",col=1:2)

plot(fit,main="応用数学",set.pars = NULL)
plot(fit22,main="遺伝統計学",set.pars = NULL)
#plot(fitall)

matplot(cbind("応用数学"=fit[[1]][,1],"統計遺伝学"=fit22[[1]][,1]),type="l",)




# additional plots
monthplot(fit)
monthplot(fit22)
monthplot(fitall)

library(forecast)
seasonplot(fit)
seasonplot(fit22)
seasonplot(fitall)


plot(blog$V1,blog$V2,type="l")

ryamada22.year <- ryamada.year <- rep(0,12)
for(i in 1:12){
	tmp <- (1:12) + (i-1) * 12
	points(blog$V1[tmp],blog$V2[tmp],col=i,pch=20)
	ryamada22.year[i] <- sum(blog$V1[tmp])
	ryamada.year[i] <- sum(blog$V2[tmp])
}

plot(ryamada22.year,ryamada.year)