# カラーパレット（透過色）
COLS <- c(rgb(1, 0, 1, .25), # ピンク
          rgb(0, 0, 1, .25), # ラベンダー
          rgb(0, 1, 0, .25)) # グリーン
GRAY <- rgb(.5, .5, .5, .25) # グレー

1 データ1

次の2つの標本X1，X2を比較に用いる。 X2はX1の約10倍の尺度を持つ。

x1 <- c(1, 2, 3, 4, 5, 5, 6, 7, 8) 
x2 <- c(50, 50, 65, 70, 80, 10, 20, 25, 40) 
hist(x1)

hist(x2)

2 四分位数

(q1 <- quantile(x1, probs = c(1/4, 2/4, 3/4)))
## 25% 50% 75% 
##   3   5   6
(q2 <- quantile(x2, probs = c(1/4, 2/4, 3/4)))
## 25% 50% 75% 
##  25  50  65

3 四分位数同士の比較（散布図）

標本の尺度が異なっていても， 分布形状が同じであれば一直線上にプロットされる。

matplot(x = q1, y = q2, pch = 16,
        main = '散布図（四分位数同士の比較）',
        xlab = '標本1分位数',
        ylab = '標本2分位数')
grid()
abline(a = 0, b = 10, col = 'red')

4 Q-Qプロット

2つの分布形状を分位数を用いて比較するための散布図を Q-Qプロット（quantile-quantile plot）という。 類似していればプロットが直線状にならぶので視覚的に分かり易い。

5 データ2

正規乱数とｔ乱数から得られた2つの標本を用いる。

set.seed(10)

n <- 100             # 標本サイズ
xn <- rnorm(n)       # 正規乱数
xt <- rt(n, df = 10) # 自由度10のｔ乱数
summary(data.frame(xn, xt))
##        xn                xt         
##  Min.   :-2.1853   Min.   :-2.0669  
##  1st Qu.:-0.8291   1st Qu.:-0.5995  
##  Median :-0.1933   Median : 0.0516  
##  Mean   :-0.1365   Mean   : 0.1114  
##  3rd Qu.: 0.5934   3rd Qu.: 0.8164  
##  Max.   : 2.2205   Max.   : 3.7912
bins <- seq(-3, 4, 1) # 階級区切

hist(xn, breaks = bins, col = COLS[1], 
     main = 'ヒストグラム 1', 
     xlab = '観測値', 
     ylab = '度数')

hist(xt, breaks = bins, col = COLS[2],
     main = 'ヒストグラム 2', 
     xlab = '観測値', 
     ylab = '度数')

6 理論分布のn分位数の計算

理論分布の面積（確率）を\(n\)分割する値（\(n\)分位数）を求める。
分割には次の2式がよく使われる。\(p_k\)は順位\(k\)までの累積確率を表す。
いずれも\(n\)が十分大きければほぼ同じ結果になるので気にしなくて良い。
（1）\(n\)分割される。ただし，分布の両脇の面積だけは\(0.5/n\)
（2）\(n+1\)等分割される。
\[ \begin{align} k = 1, 2,\cdots,n~のとき，\\ (1)\quad p_k &= \frac{k - 0.5}{n}\\ (2)\quad p_k &= \frac{k}{n+1}\\ \end{align} \]

p <- function(k) (k - 0.5) / n
k <- 1:n
q <- qnorm(p = p(k)) # 累積確率を与え分位数を求める。

理論分布が100分割される様子を次に示す。
この例では，理論分布は標準正規分布で計算式（1）を用いた。
確率密度曲線，\(x\)軸，縦線でそれぞれ囲まれる面積（確率）が\(1/n\)となる。

x <- seq(-4, 4, 0.1)
y <- dnorm(x)
matplot(x, y, type = 'l')
abline(v = q, col = COLS[1])

7 標本分布のn分位数の計算

標本サイズがｎのｎ分位数は，単に小さい順に並べ変えるだけでよい。

xn.o <- xn[order(xn)]
xt.o <- xt[order(xt)]

8 正規Q-Qプロット

\(x\)軸を標準正規分布から計算した分位数（正規分位数），\(y\)軸を標本分位数とした散布図のことを正規Q-Qプロットという。 標本分布が正規分布とみなせるか，つまり，正規性の確認を視覚的に行うためのグラフである。
なお，Q-Qプロット自体は正規分布との比較に限定するものではなく， 標本分布 VS 標本分布，理論分布 VS 理論分布などでもよい。 正規分布と比較する場合，分かり易くするため「正規」の文字を入れている。
プロットが直線状になっていれば正規分布とみなせると判断する。

matplot(x = q, y = xn.o, type = 'n',
        xlim = c(-4, 4),
        ylim = c(-4, 4),
        main = '正規Q-Qプロット',
        xlab = '正規分位数',
        ylab = '標本分位数')
matpoints(x = q, y = xn.o, pch = 16, col = COLS[1])
matpoints(x = q, y = xt.o, pch = 16, col = COLS[2])
abline(a = 0, b = 1, col = 'red')
legend('topleft', col = COLS[1:2], pch = 16,
       legend = c('標本1（正規乱数)', '標本2（ｔ乱数）'))

9 正規Q-Qプロット（Rパッケージ利用）

carパッケージを使うと信頼区間（デフォルトα＝0.05）が描画される。 この区間外のプロットがあると，有意水準5％で正規分布と有意に異なる。 （参考）

library(car)

qqPlot(xn,
       main = 'N(0, 1)',
       xlab = '正規分位数',
       ylab = '標本分位数')

## [1] 88 54
qqPlot(xt, 
       main = 't(d.f.=10)',
       xlab = '正規分位数',
       ylab = '標本分位数')

## [1] 52 49

10 Shapiro-Wilk正規性検定

帰無仮説 (H0)：標本分布は正規分布に従う

ｐ値（p-value）が0.05以下だと， 有意水準5％で帰無仮説は棄却される（正規分布ではない）。 0.05より大きいと帰無仮説を棄却できない（正規分布でないとはいえない）。

shapiro.test(xn)
## 
##  Shapiro-Wilk normality test
## 
## data:  xn
## W = 0.98908, p-value = 0.5911
shapiro.test(xt)
## 
##  Shapiro-Wilk normality test
## 
## data:  xt
## W = 0.96269, p-value = 0.006286