以下のスクリプトは、データ系列のベイジアン変化点(bcp)を実行してプロットします。ここでは、データファイル2GB01.csvの列13〜18をx<-13:18として参照する6列のデータが含まれています。問題は、列18の出力プロットのみがプロットされることです。1つのグラフィックウィンドウパーティションに各列データセットの6つのグラフをプロットするようにスクリプトを変更する必要があります。何か助けはありますか?
データファイルは次のようになります....col1col2 .... col13 col14 col15 col16 col17 col18 11 14 17 87 90 10 20 76 23 87 87 16 09 178765。。。。。。。。。。。。。。。。。。。。。。。。
スクリプトは以下のとおりです。
ここにコードを入力してください
data2<-read.csv("2GB01.csv", header=TRUE, row.names=1)
x <- 13:18
for(i in seq(along=x)){
ond<-data2[,i]
bcp.1<-bcp(as.vector(ond))
if (require("strucchange")) {
bp <- breakpoints(ond ~ 1, h = 2)$breakpoints
rho <- rep(0, length(ond))
rho[bp] <- 1
b.num<-1 + c(0,cumsum(rho[1:(length(rho)-1)]))
bp.mean <- unlist(lapply(split(ond,b.num),mean))
bp.ri <- rep(0,length(ond))
for (i in 1:length(bp.ri)) bp.ri[i] <- bp.mean[b.num[i]]
xax<-seq(1960, 2010, length=51)
op <- par(mfrow=c(2,1),col.lab="black",col.main="black")
op2 <- par(mar=c(0,4,4,4.5),xaxt="n",yaxt="n", cex.axis=0.75,las=2)
plot(xax, ond, col="grey", pch=20, xlab="",ylab="", axes=T)
lines(xax, bcp.1$posterior.mean, lwd=2)
axis(4, yaxt="s")
mtext('Posterior mean', las=0, side=4, line=3.5)
lines(bp.ri, col="blue")
par(op2)
op3 <- par(mar=c(5,4,0,4.5), xaxt="s",yaxt="n", cex.axis=0.75,las=2)
plot(xax, bcp.1$posterior.prob, xlab="Year", ylab="Posterior probability", type="l", xlim=c(1960,2010), ylim=c(0,1),las=0)
for (i in 1:length(bp.ri)) abline(v=xax[bp[i]], col="blue")
axis(2, yaxt="s")
par(op3)
par(op)
} else {
cat("strucchange is not loaded")
}}