####################################
#  第4章 モデリング
####################################
# AICによるBox-Cox変換パラメータの選択
par(mar=c(2,2,2,1)+0.1)
data(WHARD)
boxcox(WHARD)

# AIC'の値のプロット
z <- boxcox(WHARD)
x <- seq(-1,1,length=21)
plot(x,z$aic.z,col="blue",type="b",pch=20)

# 補正前のAICの値
plot(z$aic,col="red",type="b",pch=20)


# Sunspotデータとその対数変換
data(Sunspot)
plot(Sunspot,ylim=c(0,200))
plot(Sunspot,log="y")

# Box-Cox変換
boxcox(Sunspot)
#
# AIC'の値のプロット
z <- boxcox(Sunspot)
x <- seq(1,-1,length=21)
plot(x,z$aic.z,col="blue",type="b",pch=20)
plot(x,z$aic.z,col="blue",type="b",pch=20,ylim=c(2430,2490))


# 分布形状データ
x <- c(-7.99,-4.01,-1.56,-0.99,-0.93,-0.80,-0.77,-0.71,-0.42,-0.02,0.65,0.78,0.80,1.14,1.15,1.24,1.29,2.91,4.84,6.82)
y <- rep(0,length=20)
plot(x,y,ann=F,bty="l",pch=19,yaxt="n")

# 分布の同一性
x1 <- c(0.26,-1.33,1.07,1.78,-0.16,0.03,-0.79,-1.55,1.27,0.56,-0.95,0.60,0.27,1.67,0.60,-0.42,1.87,0.65,-0.75,1.52)
x2 <- c(1.70,0.84,1.34,0.11,-0.88,-1.43,3.52,2.69,2.51,-1.83)
y1 <- rep(0,length=20)
y2 <- rep(1,length=10)
plot(x1,y1,ann=F,bty="l",pch=19,yaxt="n",ylim=c(-1,2),xlim=c(-2,2))
par(new=T)
plot(x2,y2,ann=F,bty="l",pch=19,yaxt="n",ylim=c(-1,2),xlim=c(-2,2),col=2)

####################################
#  第5章 最小二乗法
####################################
#
# 多項式回帰（データの作成とモデルの推定）
 x <- seq(0,1,length=21)
 y <- c(0.125, 0.156, 0.193, -0.032, -0.075, -0.064, 0.006, -0.135, 0.105, 0.131, 0.154, 0.114, -0.094, 0.215, 0.035, 0.327, 0.061, 0.383, 0.357, 0.605, 0.499)
plot(x,y,pch=19,ylim=c(-0.2,0.6))
 polreg(y,7)
#
# test data 2
x <- seq(0,1,length=20)
z <- c(0.854,0.786,0.706,0.763,0.772,0.693,0.805,0.739,0.760,0.764,0.810,0.791,0.798,0.841,0.882,0.879,0.863,0.934,0.971,0.985)
plot(x,z,pch=19,ylim=c(0.6,1.0))
polreg(z,9)


# maxtempデータに三角関数モデルをあてはめる
data(Temperature)   # Highest Temperature Data of Tokyo
lsqr(Temperature)
#
# AIC'の値のプロット
z <- lsqr(Temperature)
x <- seq(0,21,length=22)
plot(x,z$aic,col="blue",type="b",pch=19)
plot(x,z$aic,col="blue",type="b",pch=19,ylim=c(2430,2490))


＃maxtempデータに多項式回帰モデル（最大次数７）をあてはめる
par(mar=c(2,2,3.5,1)+0.1)
data(Temperature)   # Highest Temperature Data of Tokyo
polreg(Temperature, 13)
#
# AIC'の値のプロット
z <- polreg(Temperature,13)
x <- seq(0,12,length=14)
plot(x,z$aic,col="blue",type="b",pch=19)
plot(x,z$aic,col="blue",type="b",pch=19,ylim=c(2450,2520))

#
# Whardデータに多項式回帰モデル（最大次数14）をあてはめる
  data(WHARD)  # Wholesale hardware data
  y <- log10(WHARD)
  polreg(y, 14)