.2022 Dec;19(12):1599-1611.

doi: 10.1038/s41592-022-01640-x. Epub 2022 Oct 27.

A framework for detecting noncoding rare-variant associations of large-scale whole-genome sequencing studies

Zilin Li^#^{1 2}, Xihao Li^#³, Hufeng Zhou³, Sheila M Gaynor³, Margaret Sunitha Selvaraj^{4 5 6}, Theodore Arapoglou³, Corbin Quick³, Yaowu Liu⁷, Han Chen^{8 9}, Ryan Sun¹⁰, Rounak Dey³, Donna K Arnett¹¹, Paul L Auer¹², Lawrence F Bielak¹³, Joshua C Bis¹⁴, Thomas W Blackwell¹⁵, John Blangero¹⁶, Eric Boerwinkle^{8 17}, Donald W Bowden¹⁸, Jennifer A Brody¹⁴, Brian E Cade^{5 19 20}, Matthew P Conomos²¹, Adolfo Correa²², L Adrienne Cupples^{23 24}, Joanne E Curran¹⁶, Paul S de Vries⁸, Ravindranath Duggirala¹⁶, Nora Franceschini²⁵, Barry I Freedman²⁶, Harald H H Göring¹⁶, Xiuqing Guo²⁷, Rita R Kalyani²⁸, Charles Kooperberg²⁹, Brian G Kral²⁸, Leslie A Lange³⁰, Bridget M Lin³¹, Ani Manichaikul³², Alisa K Manning^{6 33 34}, Lisa W Martin³⁵, Rasika A Mathias²⁸, James B Meigs^{5 6 36}, Braxton D Mitchell^{37 38}, May E Montasser³⁹, Alanna C Morrison⁸, Take Naseri⁴⁰, Jeffrey R O'Connell³⁷, Nicholette D Palmer¹⁸, Patricia A Peyser¹³, Bruce M Psaty^{14 41 42}, Laura M Raffield⁴³, Susan Redline^{19 20 44}, Alexander P Reiner^{29 41}, Muagututi'a Sefuiva Reupena⁴⁵, Kenneth M Rice²¹, Stephen S Rich³², Jennifer A Smith^{13 46}, Kent D Taylor²⁷, Margaret A Taub⁴⁷, Ramachandran S Vasan^{24 48}, Daniel E Weeks⁴⁹, James G Wilson⁵⁰, Lisa R Yanek²⁸, Wei Zhao¹³; NHLBI Trans-Omics for Precision Medicine (TOPMed) Consortium; TOPMed Lipids Working Group; Jerome I Rotter²⁷, Cristen J Willer^{51 52 53}, Pradeep Natarajan^{4 5 6}, Gina M Peloso^{23 24}, Xihong Lin^{54 55 56}

Collaborators, Affiliations

PMID:36303018
PMCID: PMC10008172
DOI: 10.1038/s41592-022-01640-x

A framework for detecting noncoding rare-variant associations of large-scale whole-genome sequencing studies

Zilin Li et al. Nat Methods.2022 Dec.

Abstract

Large-scale whole-genome sequencing studies have enabled analysis of noncoding rare-variant (RV) associations with complex human diseases and traits. Variant-set analysis is a powerful approach to study RV association. However, existing methods have limited ability in analyzing the noncoding genome. We propose a computationally efficient and robust noncoding RV association detection framework, STAARpipeline, to automatically annotate a whole-genome sequencing study and perform flexible noncoding RV association analysis, including gene-centric analysis and fixed window-based and dynamic window-based non-gene-centric analysis by incorporating variant functional annotations. In gene-centric analysis, STAARpipeline uses STAAR to group noncoding variants based on functional categories of genes and incorporate multiple functional annotations. In non-gene-centric analysis, STAARpipeline uses SCANG-STAAR to incorporate dynamic window sizes and multiple functional annotations. We apply STAARpipeline to identify noncoding RV sets associated with four lipid traits in 21,015 discovery samples from the Trans-Omics for Precision Medicine (TOPMed) program and replicate several of them in an additional 9,123 TOPMed samples. We also analyze five non-lipid TOPMed traits.

PubMed Disclaimer

Conflict of interest statement

Competing interests

S.M.G. is now an employee of Regeneron Genetics Center. J.B.M. is an Academic Associate for Quest Diagnostics R&D. For B.D.M.: The Amish Research Program receives partial support from Regeneron Pharmaceuticals. M.E.M. reports grant from Regeneron Pharmaceutical unrelated to the present work. B.M.P. serves on the Steering Committee of the Yale Open Data Access Project funded by Johnson & Johnson. L.M.R. is a consultant for the TOPMed Admistrative Coordinating Center (through Westat). For S.R.: Jazz Pharma, Eli Lilly, Apnimed, unrelated to the present work. The spouse of C.J.W. works at Regeneron Pharmaceuticals. P.N. reports investigator-initiated grants from Amgen, Apple, AstraZeneca, Boston Scientific, and Novartis, personal fees from Apple, AstraZeneca, Blackstone Life Sciences, Foresite Labs, Novartis, Roche / Genentech, is a co-founder of TenSixteen Bio, is a shareholder of geneXwell and TenSixteen Bio, and spousal employment at Vertex, all unrelated to the present work. X. Lin is a consultant of AbbVie Pharmaceuticals and Verily Life Sciences. The remaining authors declare no competing interests.

Figures

**Extended Data Fig. 1|. Rare variant (MAF< 0.01) distribution in the discovery phase using TOPMed cohorts (n=21,015).**
Variant categories are defined by GENCODE VEP categories.

Extended Data Fig. 2|. Manhattan plots and Q-Q plots for unconditional gene-centric noncoding analysis and sliding window analysis of high-density lipoprotein cholesterol (HDL-C) in the discovery phase (n=21,015).
a, Manhattan plots for unconditional gene-centric noncoding analysis of protein-coding gene. The horizontal line indicates a genome-wide STAAR-OP value threshold of $3.57 \times 10^{- 7}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (20, 000 \times 7) = 3.57 \times 10^{- 7})$ . Different symbols represent the STAAR-OP value of the protein-coding gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS). Promoter_CAGE and promoter_DHS are the promoters with overlap of Cap Analysis of Gene Expression (CAGE) sites and DNase hypersensitivity (DHS) sites for a given gene, respectively. Enhancer_CAGE and enhancer_DHS are the enhancers in GeneHancer predicted regions with the overlap of CAGE sites and DHS sites for a given gene, respectively.b, Quantile-quantile plots for unconditional gene-centric noncoding analysis of protein-coding gene. Different symbols represent the STAAR-OP-value of the gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS).c, Manhattan plots for unconditional gene-centric noncoding analysis of ncRNA gene. The horizontal line indicates a genome-wide STAAR-OP value threshold of $2.50 \times 10^{- 6}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / 20, 000 = 2.50 \times 10^{- 6})$ .d, Quantile-quantile plots for unconditional gene-centric noncoding analysis of ncRNA gene.e, Manhattan plot for 2-kb sliding windows. The horizontal line indicates a genome-wideP value threshold of $1.88 \times 10^{- 8}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (2.66 \times 10^{6}) = 1.88 \times 10^{- 8})$ .f, Quantile-quantile plot for 2-kb sliding windows. In panels,a,c ande, the chromosome number are indicated by the colors of dots. In all panels, STAAR-O is a two-sided test.

Extended Data Fig. 3|. Manhattan plots and Q-Q plots for unconditional gene-centric noncoding analysis and sliding window analysis of low-density lipoprotein cholesterol (LDL-C) in the discovery phase (n=21,015).
a, Manhattan plots for unconditional gene-centric noncoding analysis of protein-coding gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $3.57 \times 10^{- 7}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (20, 000 \times 7) = 3.57 \times 10^{- 7})$ . Different symbols represent the STAAR-OP-value of the protein-coding gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS). Promoter_CAGE and promoter_DHS are the promoters with overlap of Cap Analysis of Gene Expression (CAGE) sites and DNase hypersensitivity (DHS) sites for a given gene, respectively. Enhancer_CAGE and enhancer_DHS are the enhancers in GeneHancer predicted regions with the overlap of CAGE sites and DHS sites for a given gene, respectively.b, Quantile-quantile plots for unconditional gene-centric noncoding analysis of protein-coding gene. Different symbols represent the STAAR-OP-value of the gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS).c, Manhattan plots for unconditional gene-centric noncoding analysis of ncRNA gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $2.50 \times 10^{- 6}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / 20, 000 = 2.50 \times 10^{- 6})$ .d, Quantile-quantile plots for unconditional gene-centric noncoding analysis of ncRNA gene.e, Manhattan plot for 2-kb sliding windows. The horizontal line indicates a genome-wideP-value threshold of $1.88 \times 10^{- 8}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (2.66 \times 10^{6}) = 1.88 \times 10^{- 8})$ .f, Quantile-quantile plot for 2-kb sliding windows. In panels,a,c ande, the chromosome number are indicated by the colors of dots. In all panels, STAAR-O is a two-sided test.

**Extended Data Fig. 4|. Manhattan plots and Q-Q plots for unconditional gene-centric noncoding analysis and sliding window analysis of triglycerides (TG) in the discovery phase (n=21,015).**
a, Manhattan plots for unconditional gene-centric noncoding analysis of protein-coding gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $3.57 \times 10^{- 7}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (20, 000 \times 7) = 3.57 \times 10^{- 7})$ . Different symbols represent the STAAR-OP-value of the protein-coding gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS). Promoter_CAGE and promoter_DHS are the promoters with overlap of Cap Analysis of Gene Expression (CAGE) sites and DNase hypersensitivity (DHS) sites for a given gene, respectively. Enhancer_CAGE and enhancer_DHS are the enhancers in GeneHancer predicted regions with the overlap of CAGE sites and DHS sites for a given gene, respectively.b, Quantile-quantile plots for unconditional gene-centric noncoding analysis of protein-coding gene. Different symbols represent the STAAR-OP-value of the gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS).c, Manhattan plots for unconditional gene-centric noncoding analysis of ncRNA gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $2.50 \times 10^{- 6}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / 20, 000 = 2.50 \times 10^{- 6})$ .d, Quantile-quantile plots for unconditional gene-centric noncoding analysis of ncRNA gene.e, Manhattan plot for 2-kb sliding windows. The horizontal line indicates a genome-wideP-value threshold of $1.88 \times 10^{- 8}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (2.66 \times 10^{6}) = 1.88 \times 10^{- 8})$ .f, Quantile-quantile plot for 2-kb sliding windows. In panels,a,c ande, the chromosome number are indicated by the colors of dots. In all panels, STAAR-O is a two-sided test.

**Extended Data Fig. 5|. Manhattan plots and Q-Q plots for unconditional gene-centric noncoding analysis and sliding window analysis of total cholesterol (TC) in the discovery phase (n=21,015).**
a, Manhattan plots for unconditional gene-centric noncoding analysis of protein-coding gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $3.57 \times 10^{- 7}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (20, 000 \times 7) = 3.57 \times 10^{- 7})$ . Different symbols represent the STAAR-OP-value of the protein-coding gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS). Promoter_CAGE and promoter_DHS are the promoters with overlap of Cap Analysis of Gene Expression (CAGE) sites and DNase hypersensitivity (DHS) sites for a given gene, respectively. Enhancer_CAGE and enhancer_DHS are the enhancers in GeneHancer predicted regions with the overlap of CAGE sites and DHS sites for a given gene, respectively.b, Quantile-quantile plots for unconditional gene-centric noncoding analysis of protein-coding gene. Different symbols represent the STAAR-OP-value of the gene using different functional categories (upstream, downstream, UTR, promoter_CAGE, promoter_DHS, enhancer_CAGE, enhancer_DHS).c, Manhattan plots for unconditional gene-centric noncoding analysis of ncRNA gene. The horizontal line indicates a genome-wide STAAR-OP-value threshold of $2.50 \times 10^{- 6}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / 20, 000 = 2.50 \times 10^{- 6})$ .d, Quantile-quantile plots for unconditional gene-centric noncoding analysis of ncRNA gene.e, Manhattan plot for 2-kb sliding windows. The horizontal line indicates a genome-wideP-value threshold of $1.88 \times 10^{- 8}$ . The significant threshold is defined by multiple comparisons using the Bonferroni correction $(0.05 / (2.66 \times 10^{6}) = 1.88 \times 10^{- 8})$ .f, Quantile-quantile plot for 2-kb sliding windows. In panels,a,c ande, the chromosome number are indicated by the colors of dots. In all panels, STAAR-O is a two-sided test.

**Fig. 1 |. Workflow ofSTAARpipeline.**
(a) Prepare the input data ofSTAARpipeline, including genotypes, phenotypes and covariates. (b) Annotate all variants in the genome using FAVORannotator through FAVOR database and calculate the (sparse) genetic relatedness matrix. (c) Define analysis units in the noncoding genome: eight functional categories of regulatory regions, sliding windows and dynamic windows using SCANG. (d) Obtain genome-wide significant associations and perform analytical follow-up via conditional analysis.

See this image and copyright information in PMC

Comment in

STAARpipeline: an all-in-one rare-variant tool for biobank-scale whole-genome sequencing data.
[No authors listed][No authors listed]Nat Methods. 2022 Dec;19(12):1532-1533. doi: 10.1038/s41592-022-01641-w.Nat Methods. 2022.PMID:36316564No abstract available.

Cited by

Disease-specific prioritization of non-coding GWAS variants based on chromatin accessibility.
Liang Q, Abraham A, Capra JA, Kostka D.Liang Q, et al.HGG Adv. 2024 Jul 18;5(3):100310. doi: 10.1016/j.xhgg.2024.100310. Epub 2024 May 21.HGG Adv. 2024.PMID:38773771Free PMC article.
Whole genome sequence analysis of blood lipid levels in >66,000 individuals.
Selvaraj MS, Li X, Li Z, Pampana A, Zhang DY, Park J, Aslibekyan S, Bis JC, Brody JA, Cade BE, Chuang LM, Chung RH, Curran JE, de las Fuentes L, de Vries PS, Duggirala R, Freedman BI, Graff M, Guo X, Heard-Costa N, Hidalgo B, Hwu CM, Irvin MR, Kelly TN, Kral BG, Lange L, Li X, Lisa M, Lubitz SA, Manichaikul AW, Michael P, Montasser ME, Morrison AC, Naseri T, O'Connell JR, Palmer ND, Peyser PA, Reupena MS, Smith JA, Sun X, Taylor KD, Tracy RP, Tsai MY, Wang Z, Wang Y, Bao W, Wilkins JT, Yanek LR, Zhao W, Arnett DK, Blangero J, Boerwinkle E, Bowden DW, Chen YI, Correa A, Cupples LA, Dutcher SK, Ellinor PT, Fornage M, Gabriel S, Germer S, Gibbs R, He J, Kaplan RC, Kardia SLR, Kim R, Kooperberg C, Loos RJF, Viaud-Martinez KA, Mathias RA, McGarvey ST, Mitchell BD, Nickerson D, North KE, Psaty BM, Redline S, Reiner AP, Vasan RS, Rich SS, Willer C, Rotter JI, Rader DJ, Lin X; NHLBI Trans-Omics for Precision Medicine (TOPMed) Consortium; Peloso GM, Natarajan P.Selvaraj MS, et al.Nat Commun. 2022 Oct 11;13(1):5995. doi: 10.1038/s41467-022-33510-7.Nat Commun. 2022.PMID:36220816Free PMC article.
Carrier frequencies, trends, and geographical distribution of hearing loss variants in China: The pooled analysis of 2,161,984 newborns.
Feng J, Zeng Z, Luo S, Liu X, Luo Q, Yang K, Zhang G, Liu J.Feng J, et al.Heliyon. 2024 Jan 28;10(3):e24850. doi: 10.1016/j.heliyon.2024.e24850. eCollection 2024 Feb 15.Heliyon. 2024.PMID:38322914Free PMC article.
DYNATE: Localizing rare-variant association regions via multiple testing embedded in an aggregation tree.
Li X, Pura J, Allen A, Owzar K, Lu J, Harms M, Xie J.Li X, et al.Genet Epidemiol. 2024 Feb;48(1):42-55. doi: 10.1002/gepi.22542. Epub 2023 Nov 28.Genet Epidemiol. 2024.PMID:38014869Free PMC article.
CWAS-Plus: Estimating category-wide association of rare noncoding variation from whole-genome sequencing data with cell-type-specific functional data.
Kim Y, Jeong M, Koh IG, Kim C, Lee H, Kim JH, Yurko R, Kim IB, Park J, Werling DM, Sanders SJ, An JY.Kim Y, et al.medRxiv [Preprint]. 2024 Apr 15:2024.04.15.24305828. doi: 10.1101/2024.04.15.24305828.medRxiv. 2024.Update in:Brief Bioinform. 2024 May 23;25(4):bbae323. doi: 10.1093/bib/bbae323.PMID:38699372Free PMC article.Updated.Preprint.

See all "Cited by" articles

References

1. Manolio TA et al. Finding the missing heritability of complex diseases. Nature 461, 747–753 (2009). - PMC - PubMed
1. Wainschtein P et al. Assessing the contribution of rare variants to complex trait heritability from whole-genome sequence data. Nature Genetics 54, 263–273 (2022). - PMC - PubMed
1. Hernandez RD et al. Ultrarare variants drive substantial cis heritability of human gene expression. Nature genetics 51, 1349–1355 (2019). - PMC - PubMed
1. Taliun D et al. Sequencing of 53,831 diverse genomes from the NHLBI TOPMed Program. Nature 590, 290–299 (2021). - PMC - PubMed
1. Flannick J et al. Exome sequencing of 20,791 cases of type 2 diabetes and 24,440 controls. Nature 570, 71–76 (2019). - PMC - PubMed

Methods-only references

1. Chen H et al. Efficient variant set mixed model association tests for continuous and binary traits in large-scale whole-genome sequencing studies. The American Journal of Human Genetics 104, 260–274 (2019). - PMC - PubMed
1. Gazal S et al. Linkage disequilibrium–dependent architecture of human complex traits shows action of negative selection. Nature Genetics 49, 1421–1427 (2017). - PMC - PubMed
1. Li X & Li Z xihaoli/STAARpipeline: STAARpipeline_v0.9.6 Version 0.9.6 10.5281/zenodo.6871504 (2022). - DOI
1. Li X & Li Z xihaoli/STAARpipelineSummary: STAARpipelineSummary_v0.9.6 Version 0.9.6 10.5281/zenodo.6871524 (2022). - DOI
1. Li X & Li Z xihaoli/STAARpipeline-Tutorial: v0.9.6 Version 0.9.6 10.5281/zenodo.6871408 (2022). - DOI

Publication types

Actions
Actions

MeSH terms

Actions
Actions
Actions
Actions
Actions
Actions

Movatterモバイル変換

Account

Save citation to file

Email citation

Add to Collections

Add to My Bibliography

Your saved search

Create a file for external citation management software

Your RSS Feed

Full text links

Actions

Share

A framework for detecting noncoding rare-variant associations of large-scale whole-genome sequencing studies

A framework for detecting noncoding rare-variant associations of large-scale whole-genome sequencing studies

Abstract

Conflict of interest statement

Figures

Comment in

Similar articles

Cited by

References

Methods-only references

Publication types

MeSH terms

Related information

Grants and funding

LinkOut - more resources

Full Text Sources