Vol.1
前言
前几期推文其实已经出过关于横断面网络分析中不同组网络比较检验,不过写的相对比较零散,以至于后台和评论区一直有小伙伴提问关于这块的内容。索性直接再出一期关于网络比较的推文,将横断面网络分析中的网络比较检验说的更全面一些。
当然,为了照顾不懂网络比较的小伙伴,我还是要说一些关于网络比较的背景知识:
首先,我们为什么要做网络比较?很简单,因为我们想要看看不同性别、抑或是患者组和对照组的人是否在网络结构上存在显著性差异。因为病人和正常人的差异很可能通过网络结构的信息表达出来了,那么我们从网络结构的对比中很有可能找出两组人最关键的差异信息,从而为治疗或干预、预防提供新的见解。
在了解使用网络比较的目的后,我们需要知道网络比较的对象有哪些?具体来说,我们可以对两组网络的全局信息进行比较,也可以对两组网络的局部信息进行比较。例如,我们可以检验两组网络的全局强度,也可以检验两组网络中某几个我们感兴趣的节点的强度(局部强度)。
最后,在知道这些之后,我们需要去对全局信息和局部信息进行显著性检验,只有通过了显著性检验的全局强度值,节点强度值、边权重等对象,才可以认为两组之间存在显著性差异。但遗憾的是,我们只有一批数据(一组男性的数据和一组女性的数据),我们拿这批数据用来构建了男性网络和女性网络,并从两个网络中分别提取出了各一个全局强度值(也可以是其他的指标)。我们从数值上看可能是男性更大,但没有额外的数据帮助我们做统计推断。这种情况下,我们想到了之前介绍的另一种方法——bootstrap抽样(当然也可以是置换检验)。我们可以使用bootstrap的非参数有放回的随机抽样在男性和女性群体中各创建2000个样本,然后重新构建2000个网络并生成2000个全局强度值,并计算虚拟样本中男性和女性全局强度差,然后观察男性和女性全局强度差值在所有样本中的95%置信区间是否跨过0,我们就可以进行统计推断。
此外,另一个需要注意的是问题是,如果我们需要针对局部信息进行网络比较检验时,由于一个网络中存在非常多的局部信息,例如存在非常多的节点强度值和边权重值。因此,我们可能进行了非常非常多次数的两两比较,这就可能造成Alpha膨胀现象,即一次比较不犯错误的概率是0.95,两次比较不犯错误的概率是0.95*0.95=0.9025,此时犯错误的概率已经高达10%左右了,那么n次比较犯错误的概率就膨胀到了1-0.95^n。因此,我们需要进行多重比较矫正。常用的多重比较矫正方法有很多,例如Bonferroni校正、FDR矫正等。
多重比较矫正是用于控制假阳性的利器,但如果我们控制的过于严格,可能又会造成漏报,因此,在绝大多数情况下,我们其实需要的是假阳性检出率和真阳性检出率的平衡。Bonferroni校正是一种非常简单粗暴的矫正方法,比较了多少次我们就在alpha水平上除以多少。例如,我们一共进行了50次的比较,那么此时我们设置的alpha水平不再是0.05,而是0.05/50=0.001,只有达到了0.001水平的显著我们才认为其显著。FDR相较于Bonferroni更加温和一点,有时候不至于一点阳性都发现不了。具体来说,FDR矫正是针对每一个p值进行的,通过将p值利用公式q=p*n/rank转换成q值,并最终利用q值来进行显著性判断。其中n表示多重比较的次数,rank表示p值从小到大排序中所处的秩次。
好了,在了解完上述知识之后,我们就可以使用R来进行横断面网络的网络比较检验啦!
Vol.2
如何实现
##导入需要的包
library(bootnet)
library(dplyr)
library(magrittr)
library(NetworkComparisonTest)
library(psych)
library(qgraph)
###导入并查看数据
data(package = "psych")#查看psych包中自带的数据集
df <- bfi[,1:26]#使用BFI这个数据集
head(df)#查看头几行的数据,26个变量分别是大五人格的五个维度还有性别
#将gender中的1和2分别用male和female替换
df<-within(df,{
gender[df$gender==1]<-"Male"
gender[df$gender==2]<-"Female"
})
###进行网络估计
#需要注意的是,这里使用了管道符进行操作,即先对df中进行两次筛选,保留女性或男性#被试后再删除性别列,仅使用25列的大五人格数据进行网络估计。
network_male <- estimateNetwork(df %>%
filter(gender == "Male") %>%
select(-gender),
default = "EBICglasso",
corMethod = "spearman")
network_female <- estimateNetwork(df %>%
filter(gender == "Female") %>%
select(-gender),
default = "EBICglasso",
corMethod = "spearman")
plot(network_male)
plot(network_female)
#当然,如果你希望两个图的坐标位置一致,你可以考虑使用下面的代码
plot(network_male)
p_male<-plot(network_male)
plot(network_female,layout = p_male$layout)
###检验男性和女性网络全局强度差异性和网络结构不变性
set.seed(2024) # 设置置换检验的随机种子
nct_results <- NCT(network_male, network_female,
it = 100,#为了运行快设置为100,实际做可设置为1000
progressbar = T)
nct_results
#这张图给出了网络不变性检验的结果和全局强度不变性检验的结果,不显著意味着男性网络和女性网络在网络结构和全局强度上不存在显著性差异。
#如果只看这张图我们很难了解到全局强度到底是个什么东西,那么我们构建一下全局强度,大家就看的更清楚了。
S <-
abs(
sum(
abs(c(network_male$graph)) -
abs(c(network_female$graph))
)
)/2
cat("Strength difference between the two networks is:", S)
#当我们运行上述代码时,我们发现我们计算的S值和全局强度不变性检验中的S值一模一样。S值也就是全局强度,其实就是女性网络和男性网络对应的权重矩阵的绝对值相减(这里是矩阵运算),这个算的就是两个网络每条边的显著性差异,然后求这个差异的总和并将其取绝对值除以2,这样算的就是网络中所有边权重的绝对值之和,在这里被叫做全局强度。
#同样的道理。我们也可以去了解什么是网络结构全局不变性。
M <-
max(
abs(c(network_male$graph) - c(network_female$graph))
)
cat("The biggest edge difference is:", M)
#这里计算的M和网络结构不变性检验中的M值一模一样,这就说明了网络结构不变性其实就是男性和女性网络中边权值差异的绝对值的最大值。
###检验男性和女性网络中边权重的差异性
nct_test_edges <- NCT(network_male, network_female,
it = 1000,
test.edges = T,
p.adjust.methods = "BH",
progressbar = T)
nct_test_edges
#这里列出了所有边之间两两比较之后进行多重比较矫正的P值。
#查找具有显著性差异的边并给出差异值
difference_value <- function(NCT, alpha = 0.05){
diff_edges <- NCT$einv.pvals %>% dplyr::filter(`p-value` <= alpha)
for (i in 1:nrow(diff_edges)) {
var_1 <- as.character(diff_edges[i, 1])
var_2 <- as.character(diff_edges[i, 2])
value_net_1 <- NCT$nw1[var_1, var_2]
value_net_2 <- NCT$nw2[var_1, var_2]
abs_difference <- abs(value_net_1 - value_net_2)
p_value <- diff_edges$`p-value`[i]
cat("Test Edge", i, "\n----\n")
cat(var_1, "and", var_2)
cat("\nNetwork 1:", value_net_1,
"\nNetwork 2:", value_net_2)
cat("\nAbsolute difference:", abs_difference,
"with p-value =", p_value, "\n----\n")
}
}#这一步构造了一个函数,用于筛选出P值显著的比较对,并计算边权重的差值是多少
difference_value(nct_test_edges)
#当我们运行之后发现出现这样的报错,这是什么原因呢?于是我们可以按照这个函数的逻辑来重新构建。
diff_edges <- nct_test_edges$einv.pvals %>% dplyr::filter(`p-value` <= 0.05)#提取显著边对
#在工作空间中,我们发现diff_edges观测值为0,这意味着并没有存在显著差异的边对,到这里我们就理解了上述报错了。
###节点中心性指标的网络比较
nct_test_centrality <- NCT(network_male, network_female,
it = 1000, test.centrality = T,
p.adjust.methods = "BH",
centrality = c("closeness",
"betweenness",
"strength",
"expectedInfluence"),
progressbar = T)
nct_test_centrality$diffcen.pval
#观察P值后我们发现均不显著。
参考文献:
Burger, J., Isvoranu, A. M., Lunansky, G., Haslbeck, J. M. B., Epskamp, S., Hoekstra, R. H. A., Fried, E. I., Borsboom, D., & Blanken, T. F. . Reporting standards for psychological network analyses in cross-sectional data. https://psyarxiv.com/4y9nz/
Epskamp, S., Cramer, A. O. J, Waldorp, L. J., Schmittmann, V. D., Borsboom, D. (2012). qgraph: Network visualizations of relationships in psychometric data. Journal of Statistical Software, 48(4), 1–18. https://doi.org/10.18637/jss.v048.i04
Epskamp, S., Borsboom, D., & Fried, E. I. (2017). Estimating psychological networks and their accuracy: A tutorial paper. Behavior Research Methods, 50(1), 195–212. https://doi.org/10.3758/s13428-017-0862-1
van Borkulo, C. D., Boschloo, L., Kossakowski, J. J., Tio, P., Schoevers, R. A., Borsboom, D., & Waldorp, L. J. (2017). Comparing network structures on three aspects: A permutation test. Manuscript submitted. https://doi.org/10.13140/RG.2.2.29455.38569
参考链接:https://reisrgabriel.com/blog/2021-10-01-nct/
PSYCH统计实验室
通知公告
网络分析课程目前开放视频课啦!
单次课200元/讲(学生),250元/讲(非学生)
共有四讲内容:
①横断面网络分析简介与基础
②网络分析与因子分析
③交叉滞后网络分析
④时间序列网络分析
购买后开放视频权限14天,可多次申请。
并赠送所有课程相关资料(无PPT)
如果想申请购买,请联系M18812507626
更多资讯
关注我们
文稿:摘星
排版:Peruere
责编:Wink
审核:摘星
本文由“Psych统计自习室”课题组原创,欢迎转发至朋友圈。如需转载请联系后台,征得作者同意后方可转载。