2変数のt検定

　1変数のt検定はよく知られているが、2変数以上の場合はどうすればよいのだろう。2変数の場合はF検定をすることも可能である。3変数以上場合は、それぞれのt検定することはよくされていると思うが、個人的な考えで、同時に検定しなければならない。

#線形モデル
$$y_t\={\beta_1}{x_{1,t}}+{\beta_2}{x_{2,t}}+\varepsilon_t,\varepsilon_t\sim{N(0,\sigma^2)}$$

####回帰係数の推定量
```math
{\begin{bmatrix} \hat{\beta_1} \\ \hat{\beta_2} \end{bmatrix}}
=
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{-1}
{\begin{bmatrix} \sum{{x_{1,t}}{y_t}} \\ \sum{{x_{2,t}}{y_t}} \end{bmatrix}}
```
####推定量と真の値の差

```math
{\begin{bmatrix} \hat{\beta_1} \\ \hat{\beta_2} \end{bmatrix}}
-{\begin{bmatrix} \beta_1 \\ \beta_2 \end{bmatrix}}
=
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{-1}
{\begin{bmatrix} \sum{{x_{1,t}}{\varepsilon_t}} \\ \sum{{x_{2,t}}{\varepsilon_t}} \end{bmatrix}}
```

####期待値と分散
```math
\begin{align}
E{\begin{bmatrix} \hat{\beta_1}-\beta_1 \\ \hat{\beta_2}-\beta_2 \end{bmatrix}}
&=
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{-1}
E{\begin{bmatrix} \sum{{x_{1,t}}{\varepsilon_t}} \\ \sum{{x_{2,t}}{\varepsilon_t}} \end{bmatrix}} \\
&=
{\begin{bmatrix} 0 \\ 0 \end{bmatrix}}
\end{align}
```
```math
\begin{align}
var{\begin{bmatrix} \hat{\beta_1}-\beta_1 \\ \hat{\beta_2}-\beta_2 \end{bmatrix}}
&=
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{-1}
var{\begin{bmatrix} \sum{{x_{1,t}}{\varepsilon_t}} \\ \sum{{x_{2,t}}{\varepsilon_t}} \end{bmatrix}} \\
&=
{\sigma}^2
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{-1}
\end{align}
```
####分散既知の場合なら、z検定
```math
{\sigma}^{-1}
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{\frac{1}{2}}
\Bigl(
{\begin{bmatrix} \hat{\beta_1}-\beta_1 \\ \hat{\beta_2}-\beta_2 \end{bmatrix}}
\Bigl)
\sim
{N\Bigl(0 , I \Bigl)}
```

####分散未知の場合なら、t検定
```math
{\hat{\sigma}}^{-1}
{\begin{bmatrix} \sum{(x_{1,t})^2} & \sum{{x_{1,t}}{x_{2,t}}} \\ \sum{{x_{1,t}}{x_{2,t}}} & \sum{(x_{2,t})^2} \end{bmatrix}}^{\frac{1}{2}}
\Bigl(
{\begin{bmatrix} \hat{\beta_1}-\beta_1 \\ \hat{\beta_2}-\beta_2 \end{bmatrix}}
\Bigl)
\sim
{t_{n-2}}
```

#シミュレーション
モンテカルロ実験
実行回数：10000回
### サンプルサイズ　n=30
|t|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|平均|-0.0090|-0.0004|

|t共分散|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|$\hat{\beta_1}$|1.2086|-0.0293|
|$\hat{\beta_2}$|-0.0293|1.0825|

|z|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|平均|-0.0077|0.0019|

|z共分散|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|$\hat{\beta_1}$|1.0773|-0.0284|
|$\hat{\beta_2}$|-0.0284|0.9622|

|有意水準|t値|$\hat{\beta_1}$棄却率|$\hat{\beta_2}$棄却率|z値|$\hat{\beta_1}$棄却率|$\hat{\beta_2}$棄却率|
|:--:|:--:|:--:|:--:|:--:|:--:|:--:|
|0.90|0.1268|0.9025|0.8985|0.1257|0.9005|0.8965|
|0.80|0.2558|0.8070|0.7991|0.2533|0.8018|0.7980|
|0.70|0.3893|0.7077|0.6972|0.3853|0.7029|0.6922|
|0.60|0.5304|0.6136|0.5990|0.5244|0.6066|0.5932|
|0.50|0.6834|0.5193|0.5015|0.6745|0.5128|0.4934|
|0.40|0.8546|0.4210|0.4020|0.8416|0.4154|0.3932|
|0.30|1.0560|0.3265|0.3035|1.0364|0.3159|0.2942|
|0.20|1.3125|0.2299|0.2070|1.2816|0.2197|0.1971|
|0.10|1.7011|0.1241|0.1033|1.6449|0.1148|0.0948|
|0.05|2.0484|0.0657|0.0504|1.9600|0.0600|0.0450|
|0.02|2.4671|0.0260|0.0192|2.3263|0.0252|0.0156|
|0.01|2.7632|0.0127|0.0078|2.5758|0.0125|0.0081|

### サンプルサイズ　n=300
|t|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|平均|-0.0023|-0.0075|

|t共分散|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|$\hat{\beta_1}$|1.0196|-0.0184|
|$\hat{\beta_2}$|-0.0184|1.0059|

|z|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|平均|-0.0027|-0.0077|

|z共分散|$\hat{\beta_1}$|$\hat{\beta_2}$|
|:--:|:--:|:--:|
|$\hat{\beta_1}$|1.0070|-0.0176|
|$\hat{\beta_2}$|-0.0176|0.9937|

|有意水準|t値|$\hat{\beta_1}$棄却率|$\hat{\beta_2}$棄却率|z値|$\hat{\beta_1}$棄却率|$\hat{\beta_2}$棄却率|
|:--:|:--:|:--:|:--:|:--:|:--:|:--:|
|0.90|0.1258|0.9007|0.9014|0.1257|0.9004|0.9004|
|0.80|0.2536|0.8039|0.8013|0.2533|0.8045|0.8009|
|0.70|0.3857|0.7016|0.6957|0.3853|0.7009|0.6943|
|0.60|0.5250|0.5962|0.5950|0.5244|0.5963|0.5957|
|0.50|0.6753|0.4999|0.4970|0.6745|0.4997|0.4961|
|0.40|0.8428|0.4019|0.3947|0.8416|0.4005|0.3942|
|0.30|1.0382|0.3006|0.2960|1.0364|0.3000|0.2967|
|0.20|1.2844|0.2020|0.2026|1.2816|0.2005|0.2008|
|0.10|1.6500|0.1011|0.1025|1.6449|0.0991|0.1011|
|0.05|1.9680|0.0532|0.0484|1.9600|0.0514|0.0488|
|0.02|2.3389|0.0206|0.0198|2.3263|0.0207|0.0192|
|0.01|2.5924|0.0107|0.0103|2.5758|0.0104|0.0102|

#結論
○　2変数のt検定を同時にしても、シミュレーション上の結果で1変数のt検定と性質が変わらない。
○　単位根を持つ過程は別の検定法にご参考ください。

##コード
```r
t <- 300 # データの数
n <- 10000 # 実行回数
b1 <- 0.4
b2 <- 0.2
b <- as.matrix(c(0.4,0.2))
X1 <- matrix(rexp(t*n,0.01),t,n)
X2 <- matrix(runif(t*n,0,100),t,n)
error <- matrix(rnorm(t*n,0,1),t,n)
Y <- b1*X1+b2*X2+error
Bvalue <- matrix(0,n,4)
for (i in 1:n) {
  y <- Y[,i]-mean(Y[,i])
  x <- cbind( X1[,i]-mean(X1[,i]), X2[,i]-mean(X2[,i]) )
  bhat <- solve(t(x) %*% x) %*% t(x) %*% y
  shat <- (t(y-x%*%bhat)%*%(y-x%*%bhat)/(t-2))^0.5
  xsq <- t(x)%*%x
  xroot <- eigen(xsq)$vectors %*% (diag(eigen(xsq)$values))^0.5
  Bvalue[i,1:2] <- shat^(-1)*c( xroot %*% (bhat-b) )
  Bvalue[i,3:4] <- c( xroot %*% (bhat-b) )
}
Average <- apply(Bvalue,2,mean)
tCovar <- cbind( c( var(Bvalue[,1]),cov(Bvalue[,1],Bvalue[,2]) ),
                 c( cov(Bvalue[,1],Bvalue[,2]),var(Bvalue[,2])) )
zCovar <- cbind( c( var(Bvalue[,3]),cov(Bvalue[,3],Bvalue[,4]) ),
                 c( cov(Bvalue[,3],Bvalue[,4]),var(Bvalue[,4])) )
Alpha <- matrix(0,12,7)
colnames(Alpha) <- c("有意水準","t","b1のt棄却率","b2のt棄却率",
                  "z","b1のz棄却率","b2のz棄却率")
Alpha[,1] <- c(0.9,0.8,0.7,0.6,0.5,0.4,0.3,
               0.2,0.1,0.05,0.02,0.01)
Alpha[,2] <- qt(1-Alpha[,1]/2,t-2)
Alpha[,5] <- qnorm(1-Alpha[,1]/2)
for(j in 1:12){
  Alpha[j,3] <- length(which(abs(Bvalue[,1])>Alpha[j,2]))/n
  Alpha[j,4] <- length(which(abs(Bvalue[,2])>Alpha[j,2]))/n
  Alpha[j,6] <- length(which(abs(Bvalue[,3])>Alpha[j,5]))/n
  Alpha[j,7] <- length(which(abs(Bvalue[,4])>Alpha[j,5]))/n
}
Alpha
```