R 中的 Marascuilo 程序
Marascuilo procedure in R
我正在进行 marascuilio 程序以比较比例之间的差异。
我正在使用以下代码(从 this tutorial:
复制和改编
## Set the proportions of interest.
p = c(0.3481, 0.1730, 0.4788)
N = length(p)
value = critical.range = c()
## Compute critical values.
for (i in 1:(N-1))
{ for (j in (i+1):N)
{
value = c(value,(abs(p[i]-p[j])))
critical.range = c(critical.range,
sqrt(qchisq(.95,3))*sqrt(p[i]*(1-p[i])/12000 + p[j]*(1-p[j])/12000))
}
}
round(cbind(value,critical.range),3)
我需要输出也打印类别的标签(例如,哪些类别正在被比较)。
因此,如果类别在单独的向量中列出,例如categories <- c("cat1", "cat2", cat"3)
,比较是cat1-cat2
,cat1-cat3
,和cat2-cat3
。
如何将这些标签附加到我的输出中?
value critical.range
[1,] 0.175 0.016
[2,] 0.131 0.018
[3,] 0.306 0.016
试试这个:
## Set the proportions of interest.
p = c(0.3481, 0.1730, 0.4788)
N = length(p)
value = critical.range = tag = c()
categories <- c("cat1", "cat2", "cat3")
## Compute critical values.
for (i in 1:(N-1)){
for (j in (i+1):N){
value <- c(value,(abs(p[i]-p[j])))
critical.range = c(critical.range,
sqrt(qchisq(.95,N-1))*sqrt(p[i]*(1-p[i])/12000 + p[j]*(1-p[j])/12000))
tag = c(tag, paste(categories[i], categories[j], sep = "-"))
}
}
df <- as.data.frame(cbind(value,critical.range, tag), stringsAsFactors = F)
df$value <- round(as.numeric(df$value),3)
df$critical.range <- round(as.numeric(df$critical.range),3)
输出:
value critical.range tag
1 0.175 0.016 cat1-cat2
2 0.131 0.018 cat1-cat3
3 0.306 0.016 cat2-cat3
请注意计算 critical.range (12000) 时的分母...它基于每个类别的样本量 - 如果每个类别没有 12000 个观测值,则需要调整 - 如果您的观察值远少于 12000,则您的临界值可能远低于该函数给您的值(因此,您应该有更少的符号。差异)。
我正在进行 marascuilio 程序以比较比例之间的差异。 我正在使用以下代码(从 this tutorial:
复制和改编## Set the proportions of interest.
p = c(0.3481, 0.1730, 0.4788)
N = length(p)
value = critical.range = c()
## Compute critical values.
for (i in 1:(N-1))
{ for (j in (i+1):N)
{
value = c(value,(abs(p[i]-p[j])))
critical.range = c(critical.range,
sqrt(qchisq(.95,3))*sqrt(p[i]*(1-p[i])/12000 + p[j]*(1-p[j])/12000))
}
}
round(cbind(value,critical.range),3)
我需要输出也打印类别的标签(例如,哪些类别正在被比较)。
因此,如果类别在单独的向量中列出,例如categories <- c("cat1", "cat2", cat"3)
,比较是cat1-cat2
,cat1-cat3
,和cat2-cat3
。
如何将这些标签附加到我的输出中?
value critical.range
[1,] 0.175 0.016
[2,] 0.131 0.018
[3,] 0.306 0.016
试试这个:
## Set the proportions of interest.
p = c(0.3481, 0.1730, 0.4788)
N = length(p)
value = critical.range = tag = c()
categories <- c("cat1", "cat2", "cat3")
## Compute critical values.
for (i in 1:(N-1)){
for (j in (i+1):N){
value <- c(value,(abs(p[i]-p[j])))
critical.range = c(critical.range,
sqrt(qchisq(.95,N-1))*sqrt(p[i]*(1-p[i])/12000 + p[j]*(1-p[j])/12000))
tag = c(tag, paste(categories[i], categories[j], sep = "-"))
}
}
df <- as.data.frame(cbind(value,critical.range, tag), stringsAsFactors = F)
df$value <- round(as.numeric(df$value),3)
df$critical.range <- round(as.numeric(df$critical.range),3)
输出:
value critical.range tag
1 0.175 0.016 cat1-cat2
2 0.131 0.018 cat1-cat3
3 0.306 0.016 cat2-cat3
请注意计算 critical.range (12000) 时的分母...它基于每个类别的样本量 - 如果每个类别没有 12000 个观测值,则需要调整 - 如果您的观察值远少于 12000,则您的临界值可能远低于该函数给您的值(因此,您应该有更少的符号。差异)。