Cell-type deconvolution of bulk-blood RNA-seq reveals biological insights into neuropsychiatric disorders

Abstract

Genome-wide association studies (GWASs) have uncovered susceptibility loci associated with psychiatric disorders such as bipolar disorder (BP) and schizophrenia (SCZ). However, most of these loci are in non-coding regions of the genome, and the causal mechanisms of the link between genetic variation and disease risk is unknown. Expression quantitative trait locus (eQTL) analysis of bulk tissue is a common approach used for deciphering underlying mechanisms, although this can obscure cell-type-specific signals and thus mask trait-relevant mechanisms. Although single-cell sequencing can be prohibitively expensive in large cohorts, computationally inferred cell-type proportions and cell-type gene expression estimates have the potential to overcome these problems and advance mechanistic studies. Using bulk RNA-seq from 1,730 samples derived from whole blood in a cohort ascertained from individuals with BP and SCZ, this study estimated cell-type proportions and their relation with disease status and medication. For each cell type, we found between 2,875 and 4,629 eGenes (genes with an associated eQTL), including 1,211 that are not found on the basis of bulk expression alone. We performed a colocalization test between cell-type eQTLs and various traits and identified hundreds of associations that occur between cell-type eQTLs and GWASs but that are not detected in bulk eQTLs. Finally, we investigated the effects of lithium use on the regulation of cell-type expression loci and found examples of genes that are differentially regulated according to lithium use. Our study suggests that applying computational methods to large bulk RNA-seq datasets of non-brain tissue can identify disease-relevant, cell-type-specific biology of psychiatric disorders and psychiatric medication.

Keywords: cell type; deconvolution; eQTL; gene expression; neuropsychiatric.

Graphical abstract

Figure 1

Cell-type-specific expression from computational deconvolution methods (A) Cell-type-proportion predictions from CIBERSORTx. A violin plot showing the range of estimated cell-type proportions for all 1,730 individuals in each of the eight major cell types. (B) R² of expression between each cell type. A heatmap of correlations (measured by R² of mean expression across samples) between the eight main cell types captured by CIBERSORTx.

Figure 2

eGenes per cell type and correlations between effect sizes (A) Number of eGenes with a significant association identified for the eight major cell types detected by CIBERSORTx; an FDR cutoff of 0.05 was used. (B) Comparison of effect size between shared cis associations with neutrophils. Restricting eGenes to those with a significant association in both the bulk eQTL analysis and neutrophil eQTL analysis, we compared the estimated effect sizes of the most significant eQTL associations. (C) Comparison of effect size between shared cis associations with monocytes. Restricting eGenes to those with a significant association in both the bulk eQTL analysis and the monocyte eQTL analysis, we compared the estimated effect sizes of the most significant eQTL associations. (D) Comparison, using reference single-cell RNA-seq, of effect sizes between shared cis associations. Restricting eGenes to those with a significant association in both the BLUEPRINT reference neutrophil eQTL analysis and our neutrophil eQTL analysis, we compared the estimated effect sizes of the most significant eQTL associations.

Figure 3

Colocalization and enrichment analyses of cell-type-specific eQTLs (A) (Top) Number of genes with coloc PP4 > 0.8 across contexts in neuropsychiatric traits. (Bottom) Number of genes with coloc PP4 > 0.8 across contexts in blood-based traits. (B) Conditional analysis of HTR6 expression in memory B cells. All genes in the locus are included in the top panel; marginally TWAS-associated genes are highlighted in blue, and those jointly significant (HTR6) are in green. The bottom panel includes a Manhattan plot of the GWAS data before (gray) and after (blue) conditioning on the imputed expression of HTR6 in memory B cells. Imputed expression of HTR6, including 238 cis-SNPs in a LASSO regression model, was obtained. Figure generated by FUSION.post_process.R script.

Figure 4

Lithium user vs. non-user analyses (A) Boxplots showing the normalized expression of KITLG (Ensembl: ENSG00000049130) in naïve B cells, stratified by dosage of SNP rs73207047 in lithium users versus nonusers. Median values are shown as a line in the box; whiskers of boxplots are 1.5 times the interquartile range. (B) Boxplots showing the normalized expression of TNFRSF11A (Ensembl: ENSG00000141655) in monocytes, stratified by dosage of SNP rs79143095 in lithium users versus nonusers. Median values are shown as a line in the box; whiskers of boxplots are 1.5 times the interquartile range. (C) Differential gene expression results for lithium users vs. lithium non-users. (Top) Volcano plot that highlights differentially expressed genes (FDR < 0.05) in red (n = 100 total differentially expressed genes). (Bottom) Average expression of each gene vs. the log fold-change (logFC) of each gene; differentially expressed genes are highlighted in red.