必要なライブラリ

library(ltm)
library(irtoys)
library(stringr)

関数

#引数 N:乱数の数、Mu:平均、Sd:標準偏差、R:相関係数
#返り値: N行2列の行列(それぞれの列が相関を持つ乱数ペア)
CorRand <- function(N, Mu, Sd, R){
 x <- rnorm(n=N, mean=Mu, sd=Sd)
 e1 <- rnorm(n=N, mean=Mu, sd=Sd)
 e2 <- rnorm(n=N, mean=Mu, sd=Sd)
 return(cbind(sqrt(R)*x+sqrt(1-R)*e1, sqrt(R)*x+sqrt(1-R)*e2))
}

数値を入れる

d   <- 1.57
#r <- 0.28
#d <- (2*r) / sqrt(1 - r^2)

min <- d/sqrt(2)
max <- -3.5
i   <- 200
mu  <- 0
sd  <- 1
control <- c("Control")
intervention <- c("Intervention")

効果量d

#quartz(type="pdf", file="CTE_d.pdf", width=5, height=5) 
par(family = "HiraKakuProN-W3")
par(mar = c(4, 4, 1, 1)) #下左上右の順

curve(dnorm(x, d, 1), -4, 4, bty = "n", col = "black", xlab = "", ylab = "", xaxt = "n", yaxt = "n", ylim = c(0, 0.6), xlim = c(-4.0, 4.0))

par(new=T)

curve(dnorm(x, 0, 1), -4, 4, bty = "n", col = "black", xlab = "SD", ylab = "", xaxt = "n", yaxt = "n", ylim = c(0, 0.6), xlim = c(-4.0, 4.0), lty=2)

axis(side=1, at=-4:4, labels = c("-4","-3","-2","-1","0","1","2","3","4"))

lines(c(0, 0),   c(0, dnorm(0, 0, 1) + 0.05), col = "black", lty = 2)
lines(c(d, d), c(0, dnorm(0, 0, 1) + 0.05), col = "black", lty = 1)

text(d/2, 0.48, paste(sprintf("d=%0.2f", d,"SD")), adj=0.5)
## Warning in sprintf("d=%0.2f", d, "SD"): one argument not used by format
## 'd=%0.2f'
text(d/2, 0.45, paste("→"), adj=0.5)

text(-1.50, 0.30, paste(control), adj=0.5)
text( d + 1.50, 0.30, paste(intervention), adj=0.3)

#dev.off()

優越率

#quartz(type="pdf", file="CTE_CLE.pdf", width=5, height=5) 
par(family = "HiraKakuProN-W3")
#par(oma = c(1,1,0,1))   #下左上右の順
par(mar = c(4, 4, 1, 1)) #下左上右の順

curve(dnorm(x, mu, sd), -4, 4, bty = "n", col = "black", xlab = "", ylab = "", xaxt = "n", yaxt= "n", ylim = c(0, 0.6), xlim = c(-3.5, 3.5))
axis(side=1, at=-4:4, labels = c("-4","-3","-2","-1","0","1","2","3","4"))

## Polygon作図用のデータ作成
## 関数部分
xx <- seq(min, max, length=i) #min〜maxまでをi等分したデータを作成
yy <- dnorm(xx, mu, sd) #変数にベクトルを入れると関数の結果もベクトルで出力される

## 下限・上限を付加
xx<- c(min, xx, max, min)
yy<- c(0, yy, 0, 0)

## Polygonで作図
polygon(xx, yy, col="black")

## 積分範囲を示す線
lines(c(min, min), c(0, dnorm(min, mu, sd)+0.05), col="black")

## 積分範囲を示す文字(srt=90:90度回転、adj=0:アライメント左揃え(0.5中央揃え、1:右揃え))
text(min, dnorm(min, mu, sd)+0.10, 
     expression(paste(Phi, ~bgroup("(",frac(d, sqrt(2)),")"))), adj=0.5)

text(-1.1, 0.10, paste("優越率"), adj=0, col="white")
text(-1.0, 0.05, paste(round(pnorm(d/sqrt(2))*100),"%"), adj=0, col="white")

#dev.off()

相関 散布図

# dを相関に
r = d/sqrt(d^2 + 4)
#r = r
#相関を持つ乱数ペア(r=1.00)
mt <- CorRand(1000, 0, 1, r)

# Y = 1 + 2 * X
x <- mt[,1] #Memo: 関数「CorRand」内とは別物(関数内の変数は外には影響しない)
y <- 1 + 2 * mt[,2]

#図化:乱数ペアのプロット 青色破線は回帰直線(Y = 1 + 2 * X)
#quartz(type = "pdf", file = "CTE_Scatter.pdf", width = 6, height = 6) 
par(family = "HiraKakuProN-W3")#Macintoshの場合
par(mar = c(4, 4, 1, 1))
plot(x, y, pch = 1, xlab = "X", ylab="Y", main="", 
     xlim = c(-3, 3), ylim = c(-5, 6),
     mgp = c(2, 0.7, 0), xaxt = "n", yaxt = "n")
axis(side = 1, at = -3:3, labels =F)
axis(side = 2, at = -6:8, labels =F)
legend("topleft", legend = sprintf("r=%0.2f", cor(x, y)), bty = "n")
lines(c(-5, 5), c(-5*2+1, 5*2+1), lty = 2, col = 4)
# 回帰直線（赤色）
lm.coef <- lm(y~x)
abline(lm.coef, col = 2, lwd = 2)

#dev.off()