Single-cell analysis of chromatin and expression reveals age- and sex-associated alterations in the human heart

Abstract

Sex differences and age-related changes in the human heart at the tissue, cell, and molecular level have been well-documented and many may be relevant for cardiovascular disease. However, how molecular programs within individual cell types vary across individuals by age and sex remains poorly characterized. To better understand this variation, we performed single-nucleus combinatorial indexing (sci) ATAC- and RNA-Seq in human heart samples from nine donors. We identify hundreds of differentially expressed genes by age and sex and find epigenetic signatures of variation in ATAC-Seq data in this discovery cohort. We then scale up our single-cell RNA-Seq analysis by combining our data with five recently published single nucleus RNA-Seq datasets of healthy adult hearts. We find variation such as metabolic alterations by sex and immune changes by age in differential expression tests, as well as alterations in abundance of cardiomyocytes by sex and neurons with age. In addition, we compare our adult-derived ATAC-Seq profiles to analogous fetal cell types to identify putative developmental-stage-specific regulatory factors. Finally, we train predictive models of cell-type-specific RNA expression levels utilizing ATAC-Seq profiles to link distal regulatory sequences to promoters, quantifying the predictive value of a simple TF-to-expression regulatory grammar and identifying cell-type-specific TFs. Our analysis represents the largest single-cell analysis of cardiac variation by age and sex to date and provides a resource for further study of healthy cardiac variation and transcriptional regulation at single-cell resolution.

Fig. 1. Overview of datasets.

a Experimental setup. b UMAP of RNA-Seq data, coloring be cell type. c Plot of marker gene expression by cell type, showing both mean expression (dot color) and the proportion of each cell type with nonzero expression for that marker (dot size). d UMAP of a co-embedding of snRNA-and ATAC-Seq data, coloring by data type. e UMAP of co-embedded RNA and ATAC data, coloring by the assigned cell type. f Enrichment of TF motifis in the accessible peaks of ATAC-Seq cells.

Fig. 2. Sex- and age-associated gene expression changes in the heart.

a Differentially expressed genes as a function of age (green) or sex (purple), FDR- 0.1. Testing used a mixed effect model to control for variation by anatomical site, cell UMI, and inter-donor variation, multiple-testing corrected using the Benjamini–Hochberg method. Total donors n = 8. b G5EA by biological pathways as a function of sex. Greer = pathway components up-regulated in females, orange = up-regulated in males. Total donors n = 8. Asterisk number indicates the FDR threshold of significance. c TF motifs enriched or depleted in accessible chromatin as a function of sexus-ing a mixed effect model controlling for age, sex, anatomical site, read depth, and inter-donor variation (see “Methods” section), FDR = 0.1. Total donors n = 9, d GSFA by biological pathways as a function of sex. Red = pathway up-regulated with age, blue = down-regulated with age. Unique donors n = 8. e TF motifs enriched or depleted in accessible chromatin as a function of age, using the same approach as (c).

Fig. 3. Meta-analysis of sex- and age-associated expression patterns across studies.

a Age and sex of all unique donors (one point = one donor) by publication. b UMAP embedding of single-nuclei RNA-Seq data from six studies, colored by cell type. c UMAP embedding as in (b), coloring points (cells) by data source. d Number of differentially expressed genes as a function of age (green) or sex (purple), FDR = 0.1. Testing used a mixed effect model controlling for variation by anatomical site, cell UMI, data source, and inter-donor variation (see “Methods” section). Uses n = 73 donors. e GSEA by biological pathways as a function of sex, looking at sex-dependent variation as fit using a mixed effect model (see “Methods” section). Green = pathway components up-regulated in females, orange = up-regulated in males. Total donors n = 73. Asterisk number indicates FDR threshold of significance. Benjamini–Hochberg method for multiple-testing correction. f GSEA by biological pathways as a function of sex. Red = pathway up-regulated with age, blue = down-regulated with age. Unique donors n = 73. g Plot of cell-type proportions per sample by sex. Note that one point represents one technical sample, n = 73 unique donors analyzed. h Plot of the adjusted proportion of neurons in a given sample by age, colored by data source. Note that n = 73 unique donors, while in this plot one point = one biological sample (repeated samples from distinct sites). Trend line represents the predicted proportion of neurons from a beta-binomial mixed effect model (see “Methods” section) with the uncertainty range representing predictions using an age coefficient ± 2 standard errors around the fit value for the age effect. Proportions were adjusted to remove effects by sex, anatomical site, and data source based on coefficients fit for those variables (see “Methods” section). i Coefficients genes as a function of age in a mixed effect model fit in neurons. Error bars represent coefficient estimate ± 2 standard errors. Q-values listed are calculated from p-values using the Benjamini–Hochberg method. j Coefficients for genes as a function of sex (positive coefficient = higher expression in males, negative = higher expression in females) in a mixed effect model fit in cardiomyocytes. Error bars represent coefficient estimate ± 2 standard errors as in (i).

Fig. 4. Enrichment of TF motifs in the accessible peaks of fetal or adult sn ATAC-Seq.

a Enrichments in cardiomyocytes. b Enrichments in cardiac neurons. c Enrichments in vascular endothelial cells. Enrichments are shown for TFs that were statistically enriched in one or both of fetal or adult analyses (FDR = 0.1). d Enrichments in macrophages (adult) versus myeloid lineage cells (fetal).

Fig. 5. Predictive models of RNA expression.

a R² for cell-type-specific models (one point = one cell type) using varied p-value cutoffs for calling motif presence in FIMO (color) and upstream/downstream regions (in bp) with respect to geneTSS. b Test set R² for different cell types, R² is listed for models that were trained using motifs found near geneTSS (in purple) or in either TSS or linked distal sites (in green). c Model accuracy on test data for models using promoter and distal motifs (x-axis) and the proportion of that cell type in the snRNA-Seq data (y-axis). d Magnitude and direction of motif coefficients in final models fit using distal and proximal sequence.