使用C++快速计算网络抗毁性

fastnc2

  根据文献报道，从网络图中提取邻阶矩阵（0、1矩阵，0代表指标间没有边，1代表指标间有边），然后依次删除节点并求剩余矩阵的“特征值”。通过将计算结果与删除节点比例做线性回归并提取slope值，比较不同网络图之间的slope数值大小（简单说，就是看回归线下降的趋势）以反映网络的抗毁性。在实践中，一个网络图有N个节点，那么当删除M个节点时，M个节点有多个情况的节点组合，目前一些文献是通过degree、BC拓扑特征对节点排序，然后依次删除。在这里，我们用另外的做法，即随机抽取M个节点并删除，正如前面所言，M个节点有多种情况组合，因此我们提供了一个参数n，用来多次迭代删除，如我们在删除123个节点，那么我们从总节点中随机删除123个节点并求剩余节点的邻接矩阵的自然连通度，这一步骤重复n次(如1000)，最后我们将n次的结算结果取均值，这样得到的结果更具代表性。
  在实际数据分析过程中，用随机删除M个节点并迭代n次的做法具有很大的计算量。我们基于C++软件、mkl库来编写代码执行上述过程；同时， 我们将软件命名为fastnc2。
1、通过conda安装软件：

conda create -n fastnc2 -c wqssf102 fastnc2

2、或直接下载编译好的版本，解压即可用:
点击下载
解压：tar xzvf fastnc2.tar.gz
设置权限：chmod a+x fastnc2

软件的参数如下：

fastnc2 --help
 Program: FastNC2 (use c++ to calculate the natural connectivity).
  Contact: Qiusheng WU (565715597@qq.com):
  -h [ --help ]               help message
  -i [ --input ] arg          input.txt,The result from get.adjacency function
                              of igraph package
  -o [ --output ] arg         output.txt
  -s [ --step ] arg (=1)      number of threads (default 1)
  -t [ --thread ] arg (=0.80) the threshold for deletion of node (default: 0.8)
  -n [ --number ] arg (=1000)  number of iterations (default: 1000)
  -j [ --jobs ] arg (=4)       number of jobs (default: 4)

  参数进一步说明：

-i为邻接矩阵，必须给出

-o 为输出文件，必须给出

-t 为删除物种（节点的阈值），默认为0.8，即从1%删到80%

-s 为步数，默认为1，假如step为3，那么就是依次删除1个节点，4个节点，7个节点，......,即第N次是删除（N+step）个节点（N>1）
因此，当step大于1时，计算时间将大大缩短，适用于网络图节点多的情况


-n 为迭代次数，默认1000，因为我们不知道环境变化会使哪些物种消失，因此这里采用了随机抽样，建议至少迭代1000次。

-j 线程数，默认为4。

  邻接矩阵如何获得，方法之一是通过R软件的igraph包导出，如：

library(igraph)
##gg为igraph对象的网络图数据
adjtab <- as.matrix(get.adjacency(gg))
####导出01矩阵，然后使用服务器计算
write.table(adjtab,"mydir/adj_tab.txt",quote=F,row.names=F,col.names=F,sep="\t")

跑分析：

conda activate fastnc2
fastnc2 -i adj_tab.txt -j 28 -n 1000 -o ncres.txt
conda deactivate

  最终结果格式如下：

nc_index	del_number
0.579490466952	0.002941176471
0.575350980163	0.005882352941
0.571343702435	0.008823529412
0.568237121284	0.011764705882
0.565799163938	0.014705882353
0.563750376701	0.017647058824
0.562035103083	0.020588235294
0.560567979574	0.023529411765
0.559203409255	0.026470588235
0.557888396978	0.029411764706
0.556521071494	0.032352941176
0.555493312657	0.035294117647
0.554323430896	0.038235294118

nc_index为删除M个节点后剩余矩阵的“特征值”，del_number为删除节点所占总节点的百分比。

  最后用R绘图：

library(ggpmisc)
grp1 <- read.csv("ck.txt",header = T,sep="\t")
grp1$grp <- "CK"
grp2 <- read.csv("NP.txt",header = T,sep="\t")
grp2$grp <- "NP"

####将每组合并
grpnc <- rbind(grp1,grp2)
##指定图例的顺序
grpnc$grp <- factor(grpnc$grp,levels = c("CK","NP"))
ggplot(grpnc, aes(del_number, nc_index,color=grp)) +##grp为将多个网络图合并时的分组
 # geom_point() +
 geom_smooth(formula = y~x,se = FALSE,method = "lm",show.legend = F)+
 stat_poly_eq(formula = y~x,size=5,family="serif",method = "lm",
              output.type = "numeric",
              parse = T,label.x = 0.3,hjust=0,
              mapping =aes(label = paste(paste(sprintf("italic(R)^2~`=`~%.2g",after_stat(`r.squared`)),"~",
                                               ifelse(after_stat(`p.value`)<=0.001,"'***'",
                                                      ifelse(after_stat(`p.value`)<=0.01,"'**'",
                                                             ifelse(after_stat(`p.value`)<=0.05,"'*'",NA)))),
                                         sprintf( "italic(Slope)~`=`~%.3g",after_stat(b_1)),sep = "*\", \"*")))+
 labs(x="Proportion of removed nodes",
      y="Natural connectivity",
      color="Group")+
 scale_color_manual(values = c("CK"="#7570B3","NP"="#1B9E77"))+
 theme_bw()+theme(panel.grid = element_blank())+
 theme(axis.title = element_text(size = 20,colour = "black"))+
 theme(axis.text = element_text(colour = "black",size = 20,angle = 0,vjust = 0.5))+
 theme(legend.title =  element_text(size = 20,colour = "black"), 
       legend.text = element_text(size = 20,colour = "black"),
       strip.text =  element_text(size = 20,colour = "black"),
       text = element_text(family = "serif"))

关于计算速度和时间

   1、若用户的网络包含的节点比较多（如5千个节点），假如用户想要迭代1000次，当服务器的运行内存不大时，我们建议分为4个任务提交到服务器中，每个任务执行250次的迭代，最后取4个任务计算结果的均值即可，这样会大大缩减用户等待的时间。特别地，可以通过设置step参数（如step设置为3）来扩大删除节点的步长，这样耗时将显著减少。总之，step是一个很有用的参数。 我们的软件还在开发中，若你对软件的功能有需求或发现软件问题，请联系作者：565715597@qq.com
  fastnc软件目前只在Linux系统下测试，其他系统没测试过（Win系统肯定是不行）。
参考文献：

[1] Wu MH, Chen SY, Chen JW, Xue K, Chen SL, Wang XM, Chen T, Kang SC, Rui JP, Thies JE, Bardgett RD, Wang YF. Reduced microbial stability in the active layer is associated with carbon loss under alpine permafrost degradation. Proc Natl Acad Sci U S A. 2021 Jun 22;118(25):e2025321118. doi: 10.1073/pnas.2025321118.

[2] G.-s. Peng, J. Wu, Optimal network topology for structural robustness based on natural connectivity. Physica A 443, 212–220 (2016). doi: 10.1016/j.physa.2015.09.023