| Background:In the era of current omics,various high-throughput omics technologies are flourishing and becoming increasingly mature.It has been basically realized to map omics molecular markers that affecting occurrence and development of diseases to biological molecular networks,forming an integrated research framework of systems biology and network medicine.From the perspective of network medicine,complex diseases are seldom resulted from the abnormality of a single gene but rather by the combined effects of multiple genes,which often form an intricate network of interactions that control the occurrence,development,and outcome of diseases.Identifying gene co-expression networks related to complex diseases can provide a more comprehensive understanding of the underlying genetic mechanisms of diseases.However,exploring the association between gene networks and complex diseases requires obtaining both gene expression profile data and disease phenotype data from the same sample.Currently,there are many international public datasets that study the relationship between transcriptomes and diseases,but their small sample sizes often lead to low statistical power.The integration of cross-omics data,represented by transcriptome-wide association studies(TWAS),provides a framework for obtaining large-scale gene expression data.However,most existing TWAS methods are limited to single-gene quantitative traits,and the only two multi-gene TWAS methods(FOCUS and FOGS)ignore the complex network relationships among multiple genes.In addition,categorical phenotypes are common in biomedical research,and directly inputting them as quantitative traits into existing TWAS frameworks will lead to loss of information inevitably.Therefore,there is an urgent need to develop TWAS network regression statistical methods with categorical phenotypes to better provide statistical support for explaining the genetic network biology mechanisms behind complex traits.Objectives:This study combines the identification of potential biological network markers for complex diseases with cross-omics TWAS methods.For categorical phenotypes,following the modeling concept of "structural decomposition→mathematical integration",a network regression proportional odds logistic model based on a two-stage TWAS framework is proposed to detect the association between specific gene networks and categorical phenotypes,identify significant genes and gene-gene interactions related to diseases,and provide statistical new ideas and methods for exploring the underlying network mechanisms of complex diseases.Methods:A network regression model called Proportional Odds LOgistic model for NEtwork regression(PoLoNet)was developed for detecting the association between gene networks and categorical phenotypes in TWAS.PoLoNet relied on two-stage TWAS framework.First,a distribution-robust nonparametric Dirichlet process regression model(DPR)was used in eQTL study to obtain estimated effects of genotypes on gene expression,which were then used in GWAS study to obtain the predicted gene expression values and represent them as nodes in the network.Next,PoLoNet used pointwise mutual information(PMI)to characterize the complex relationships,i.e.,edges,between network nodes of multiple types.To highlight the advantages of PoLoNet in detecting complex relationships between nodes,this study compared it with the method of representing edges based on product moment(PM),which had good performance in capturing linear relationships.Finally,all nodes and edges were included in the model for association analysis.The data used in this study were all from publicly available databases,where gene expression data come from GEUVADIS study,disease network structure data come from KEGG(Kyoto Encyclopedia of Genes and Genomes),GWAS and phenotype data come from UK Biobank.1.In the statistical simulations,we used a branch of the Alzheimer’s disease pathway(hsa05010-nt06412)in KEGG for simulation studies.The AD branch network includes 12 nodes and 13 edges.Based on GEUVADIS dataset,we obtained the predicted gene expression values using both the DPR and the Bayesian sparse linear mixed model(BSLMM),then used PMI and the PM to characterize the general relationships between nodes,with the PM-based method called PPNT(Proportional odds logistic model using the Product moment in Network TWAS).We then comprehensively evaluated the type I error rate and power of PoLoNet under various simulation scenarios,including different outcome categorical types,sample sizes,sample proportions between categories(including extremely unbalanced samples),different correlation patterns between nodes,and different gene expression prediction methods.In addition,to ensure the completeness of the simulation experiment,we also included simulations of randomly selected nodes and edges as well as missing gene nodes in the network,and compared PoLoNet with single-gene TWAS method.2.In the real data analysis,blood pressure and bipolar and major depression status were chosen as the two traits from the UK Biobank,and the data types were categorized into binary and ordinal categorical phenotypes according to the guidelines.Biological networks related to these two traits were selected from the KEGG database,including 22 networks related to blood pressure and 9 networks related to bipolar and major depression status.Genes within each network were matched with eQTL data from GEUVADIS to obtain the final gene network for inclusion in the model.Based on the TWAS framework,gene expression prediction values(network nodes)were obtained.and PMI and PM were calculated to represent the network edges.Finally,the outcomes of two categories of disease were analyzed to demonstrate the advantages of PoLoNet in discovering significant disease-related genes and gene-gene interactions in real data analysis.Considering that nodes and edges in TWAS network regression were often highly correlated and the commonly used multiple testing method Bonferroni correction is too stringent,this study used false discovery rate(FDR)to adjust Pvalues and set the FDR significance threshold to 0.05.Results:1.In the various simulated scenarios with binary outcomes,both PoLoNet and PPNT had stable type I error control rates when detecting effecting nodes and effecting edges.As expected,when only evaluating node effects,the power of both methods increased as increasing sample size and as the case-control ratio became more balanced,but had little difference.When evaluating the power of edge effects,as sample size increased and the casecontrol ratio became balanced,PoLoNet had significantly higher power than PPNT with nonlinear relationships between nodes and grow faster than PPNT,even under extremely unbalanced sample ratios.In addition,the type of nonlinear relationship between nodes may also affect the power of PoLoNet specific nonlinear relationship(for example:y=sin(x)2)will make the PoLoNet had higher power.When the relationships between nodes were linear,the power of PoLoNet was slightly lower than that of PPNT,as the product moment was the gold standard for measuring linear relationships in this case.The same conclusion was reached when using BSLMM to predict gene expression.When the outcome was an ordinal categorical phenotype,predicting gene expression with DPR and BSLMM,PoLoNet and PPNT were equally stabilize controlling type I error rates and had similar power in all simulated scenarios.2.In the analysis of real data for a network of 22 genes related to blood pressure,the results of detecting effecting nodes and effecting edges by PoLoNet and PPNT were consistent with those of the simulation.When blood pressure was treated as a binary phenotype,PoLoNet successfully identified 39 genes and 71 edges,while PPNT identified 39 genes and 66 edges,with 37 overlapping genes and 63 overlapping edges.When blood pressure was treated as an ordinal categorical phenotype,PoLoNet identified 63 genes and 78 edges,while PPNT identified 58 genes and 68 edges,with 56 overlapping genes and 67 overlapping edges.Real data analysis of 9 gene networks for bipolar disorder and major depression,PoLoNet and PPNT also had comparable performance in detecting effecting nodes,but PoLoNet identified more edges than PPNT.which was consistent with the results of blood pressure.To some extent,the advantages of PoLoNet in detecting significant genes and gene-gene interactions in complex disease gene networks are explained.Conclusion:In this study,we incorporated gene networks into a two-stage TWAS framework and developed a novel network regression method called PoLoNet for detecting associations between specific networks and the target binary or ordinal categorical phenotype of interests.PoLoNet relied on DPR to obtain optimal weights for gene expression prediction and introduced PMI to capture the general correlation between nodes.More importantly,PoLoNet can simultaneously identify specific node and edge effects in the network.Nonparametric kernel density estimation was employed in PoLoNet to estimate PMI,thereby avoiding the risk of misspecification of the joint distribution of two genes in the network.Extensive simulation results demonstrated that PoLoNet captured general relationships between different genes effectively and had higher power than other methods.Furthermore,PoLoNet remained robust to different gene expression prediction models and different proportions of categorical phenotypes in TWAS.In addition,PoLoNet has high computational efficiency. |
