Transcriptome analysis uncovers the autophagy‐mediated regulatory patterns of the immune microenvironment in dilated cardiomyopathy

Abstract The relationship between autophagy and immunity has been well studied. However, little is known about the role of autophagy in the immune microenvironment during the progression of dilated cardiomyopathy (DCM). Therefore, this study aims to uncover the effect of autophagy on the immune microenvironment in the context of DCM. By investigating the autophagy gene expression differences between healthy donors and DCM samples, 23 dysregulated autophagy genes were identified. Using a series of bioinformatics methods, 13 DCM‐related autophagy genes were screened and used to construct a risk prediction model, which can well distinguish DCM and healthy samples. Then, the connections between autophagy and immune responses including infiltrated immunocytes, immune reaction gene‐sets and human leukocyte antigen (HLA) genes were systematically evaluated. In addition, two autophagy‐mediated expression patterns in DCM were determined via the unsupervised consensus clustering analysis, and the immune characteristics of different patterns were revealed. In conclusion, our study revealed the strong effect of autophagy on the DCM immune microenvironment and provided new insights to understand the pathogenesis and treatment of DCM.

including genetic mutations, infections, inflammations, autoimmune diseases, exposure to toxins and endocrine or neuromuscular causes. 1 Due to different causes, DCM is usually classified as familial (genetic) or nonfamilial (nongenetic) forms. Around 15-30% DCM patients were diagnosed as familial DCM. 3,4 In the past 20 years, with the help of next-generation sequencing, several mutations have been identified as DCM marker genes, [5][6][7] which are LMNA44, MYH7, TNNT2, TTN, RBM20, BAG3 and so on. [8][9][10][11][12] Among them, TTN truncation mutation is the most common cause of DCM, occurring in ~25% of familial DCM cases and 18% of sporadic cases. 7,13 Although genetic prediction based on gene mutation is becoming a useful tool in clinic, there are limited reports on DCM diagnosis based on gene expression.
Besides genetic mutations, immune microenvironment changes are also leading causes for DCM. 7,14 Myocardial damage caused by DCM leads to inflammation with recruited immune cells into heart to repair damaged myocardium. Pathological examination of myocardial biopsy samples (or autopsy) from DCM patients often reveals evidence of activated inflammatory cell infiltration with gene expression patterns compatible with activated immune cells. 1,15 Moreover, cardiac-specific autoantibodies can be detected from ~60% DCM patients and their relatives, which directly affect cardiomyocyte function and disease prognosis. [16][17][18] In patients, inflammatory DCM may be familial and related to HLA antigens. 19 For instance, cardiac infiltrated immune cells with abnormal class II HLA expression occupied approximately 50% among all biopsy samples. 20 Therefore, determining the proportion and type of infiltrated immune cells and HLA genes expression pattern is crucial to establish the best treatment plan.
Among the entire immune system, autophagy is essential for cell development, function and homeostasis. 21,22 Cell-autonomous inflammation is one of the key contributions of autophagy to the immunity. 23,24 Furthermore, it is also proved that autophagic activation participates in innate immunity response by mediating foreign pathogens clearance. 25 By affecting cell metabolism, cytoplasmic quality and tissue homeostasis, autophagy has various connections with many human diseases. [26][27][28] Recent years, autophagy has emerged as a major regulator of cardiac homeostasis. Autophagy preserves cardiac structure and function in normal conditions and is activated in a stress response, which will limit damage in both physiological and pathological conditions. 29,30 There are some studies focusing on the effects of autophagy on DCM. For example, Gil-Cayuela et al. analysed the autophagy-related gene expression changes in DCM patients by RNA-seq. 31 They found that in DCM patients, the expression changes of NRBP2 and CALCOCO2 were related to LV dysfunction and remodelling. Furthermore, primary fibroblasts from severe autosomal recessive DCM patients with mutations in PLEKHM2 gene exhibited abnormal endosomal subcellular distribution, abnormal lysosomal localization, and impaired autophagic flux marked by RAB5, RAB7 and RAB9 respectively. 32 In addition, the combination of overactivated AKT-mTOR pathway and defective autophagy has been described in a mouse model of DCM carrying LMNA mutations. 33 This study showed that intraperitoneal injection of temsirolimus, a derivative of the mTOR inhibitor rapamycin (sirolimus), in mice could reduce cardiac dilatation and enhance cardiac function. The authors demonstrated temsirolimus treatment decreased mTORC1 signalling, increased LC3-II and decreased p62 protein level. Although all these studies have shown that there is a link between autophagy activity and DCM disease progression, this only explains the tip of the iceberg of autophagy's impact on DCM. The mechanism of how a large number of autophagy genes affect the immune microenvironment to regulate DCM remains to be explored.
Here, we systematically evaluated the effect of autophagy on DCM immune microenvironment. The risk prediction model based on the expression of 13 DCM-related autophagy genes could effectively distinguish healthy and DCM samples. Next, to explore the connection between immune microenvironment and DCM, immunocyte infiltration, immune responses and HLA status in DCM were investigated, which showed strong correlation with autophagy. To further study how autophagy regulates DCM, the unsupervised clustering on the autophagy expression profile in DCM samples was performed, and two autophagy-mediated regulation patterns in DCM were determined. These two subtypes showed distinct characteristics of immune microenvironment. Gene set variation analysis (GSVA) and the functional enrichment analysis indicated that the two autophagy-mediated patterns are mainly respectively responsible for disease and metabolism related pathways. In addition, through weighted gene co-expression network analysis (WGCNA), a co-expression module (blue module) which strongly correlated with autophagy subtypes was identified. The findings provide a comprehensive overview about how autophagy regulates DCM through the immune microenvironment.

| Data preprocess
The data used in this study were RNA sequencing data of the LV from

healthy donors and DCM samples from The Myocardial
Applied Genomics Network (MAGNet; www.med.upenn.edu/magnet). The LV free-wall tissues were harvested at the time of cardiac surgery from subjects with heart failure undergoing transplantation and from unused donor hearts with apparently normal function.
The expression matrix was reserved in the gene expression omnibus (GEO) database (https://www.ncbi.nlm.nih.gov/geo/query/ acc. cgi?acc=GSE14 1910), 34 and the clinical characteristics of involved samples in GSE141910 were listed in Table S1. In addition, the validation of autophagy model was performed in another microarray dataset GSE57338 (https://www.ncbi.nlm.nih.gov/geo/query/ acc.cgi?acc=GSE57338). 35 The gene expression was detected by Affymetrix Human Gene 1.1 ST Array microarray. Gene probes were annotated as gene symbols. Probes without matching gene symbols and matching multiple symbols were excluded. Gene expression value of duplicate gene symbol was calculated as the median value.

| Alteration analysis of autophagy between DCM and healthy samples
The differential analysis of autophagy genes between healthy and DCM samples was carried out using R package 'limma', 36 and differentially expressed autophagy genes were identified with adjusted p-value <0.01 and |log2FC| > 0.5. The protein-protein interaction network was constructed via STRING database (https://strin g-db.org/). 37 The DCM-related autophagy genes were identified by univariate logistic regression with the cut-off criteria of p value <0.0001. Then least absolute shrinkage and selection operator (LASSO) regression was performed to select the most useful biomarkers among DCM-related autophagy genes and the diagnostic model with nonzero coefficients was constructed using the R pack- Principal components analysis (PCA) was used to reduce the number of dimensions to find similarities and differences of RS between healthy and DCM samples. The ROC curve was plotted to assess the classification performance of the classifier.

| Validation of DCM-related autophagy genes using Real-Time Quantitative PCR (RT-qPCR)
Serum samples were lyzed in QIAzol Lysis Reagent and miRNeasy Serum/Plasma Kit(QIAGEN)and quantified using Nanodrop 2000 software (Nanodrop Products). For mRNA quantification, cDNA was synthesized using the PrimeScript RT reagent kit (Takara,). RT-qPCR was performed using PCR 218 reaction mixture (11203ES08,) on a

| Correlation analysis between immune microenvironment and autophagy
Immune cell infiltration was calculated using CIBERSORT algorithm based on the expression matrix of 22 types of immunocytes. 38 The relative enrichment of the activity of immune pathways was determined by single-sample gene set enrichment analysis (ssGSEA).
The immune reaction gene sets used to evaluate the activity of immune-related pathways were obtained from ImmPort (http:// www.immpo rt.org). 39 The relative fraction of immunocytes, enrichment score of immune-related pathways and expression of HLA genes between healthy and DCM samples were compared using the Wilcox test. The correlation analysis between the relative fraction of immunocytes, enrichment score of immune pathways and expression of HLA genes and autophagy genes were done by Spearman correlation analysis with R package 'corrplot'.

| Identification of autophagy expression patterns
Unsupervised clustering analysis was performed to identify distinct autophagy expression patterns according to the expression of all 201 autophagy genes. Consensus clustering was implemented to evaluate the autophagy expression pattern using the R package 'ConsensuClusterPlus'. 40 This algorithm was repeated 1000 times to ensure the stability of classification. The expression of the DCMrelated autophagy genes, immunocyte relative fraction, immune reaction scores and expression of HLA genes between the two subtypes were compared using the Wilcox test.

| Functional enrichment analysis of the two autophagy expression patterns
To screen for autophagy expression pattern related genes, differentially expressed genes between two subtypes were identified using the R package 'limma'. 36 The criterion for DEGs was |log2FC| > 0.5 and adjusted p-value <0.001. The biological characteristics of autophagy phenotype-related genes and were uncovered by GO and KEGG enrichment analysis using the R package 'clusterProfiler'. 41 To further reflect biological changes that occurred in each subtype, the gene-sets of 'h.all.v7.0.symbols' and 'c2.cp.kegg.v7.0.symbols' were downloaded from MSigDB (http://www.gsea-msigdb. org/). 42 Next, the activity of HALLMARKS and KEGG pathway was quantified using the GSVA algorithm, 43 and was compared between two subtypes using the Wilcox test. Pathways with adjusted p-value <0.01 were considered to be significantly different in activity between the two subtypes.

| Identification of autophagy expression pattern related gene modules
Prior to analyse co-expressed genes, the coefficient of variation (CV) of each gene based on their expression values were calculated and genes with low variability (CV <0.1) were filtered. Next, the R package WGCNA was used to identify modules of highly correlated genes across all samples. 44 The parameters of WGCNA used default settings, except for the power was 2, network type was signed, min-ModuleSize was 30. Eigengenes and clusters were calculated based on the correlations to quantify the co-expression similarity of entire modules, using a strict cut-off of 0.25, corresponding to correlation of 0.75. Then, kME, known as the module membership value, was calculated using SignedKME algorithm to represent the correlation between a gene and the module eigengene value. Intramodular hub genes of module blue were identified if their kME >0.95.

| The landscape of autophagy gene alterations in DCM
To explore the expression status in DCM, 222 known autophagy  Table S3). In particular, NAMPT had the largest fold change and CALCOCO2 had the most statistically significant change. To reveal the interactions between these autophagy genes, the protein-protein interaction network was constructed ( Figure 1D). Among them, HIF1A and BCL2L1 interacted with more than 10 autophagy genes.

| Autophagy is associated with immune microenvironment in DCM
To further investigate the biological connections between autophagy genes and immune microenvironment, the correlation analysis for  Table S6), suggesting a great change of immune microenvironment during DCM progression. Correlation analysis showed autophagy genes were closely related to many immunocytes ( Figure 3D, Table S7). The significant positively correlated immunocyteautophagy gene pair is CALCOCO2-eosinophils ( Figure S2A), while the most positively correlated is NAMPT-eosinophils ( Figure S2B). In contrast, the most negatively correlated pair is CX3CL1-eosinophils ( Figure S2C).  Likewise, the activity of immune-related pathways and expression levels of HLA genes were calculated and obvious changes were observed between healthy and DCM samples ( Figures 3B,C). Interestingly, nearly all HLA genes were highly  (Table S12). PCA analysis revealed that there was a remarkable difference in transcriptome between the two subtypes ( Figure 4D).
Then, the expression of the 23 significantly dysregulated autophagy genes were compared between two DCM subtypes, and only 16 of them dramatically altered ( Figure 4E,F). This indicated that these 16 autophagy genes may be involved in different mechanisms of regulating DCM, while the rest may regulate DCM indiscriminately.
To uncover the differences of immune microenvironment characteristics between the two subtypes, infiltrating immunocytes, immune response gene-sets and HLA gene expression were evaluated. As expected, the two subtypes demonstrated very distinct autophagy-mediated immune characteristic ( Figure 5, Table S13).
For example, much more eosinophils and resting memory CD4 T cells were observed in subtype-1, while higher levels of Monocytes, CD8 T cells and regulatory T cells (Tregs) were found in subtype-2 (p < 0.0001). As for immune responses, TCR signalling pathway and TGF family member receptor related genes were enriched more in subtype-1, while antigen processing and presentation and TNF family members receptors were more active in subtype-2 (p < 0.0001).
In particular, with the exception of HLA-F-AS1, all dysregulated HLA genes had higher expression levels in subtype-2, suggesting subtype-2 represents immune enrichment.

| Biological functions behind autophagy expression patterns
Differences in expression and immune microenvironment characteristics illustrated that there are two mechanisms for autophagy genes to regulate DCM. In order to clarify the respective roles of the two subtypes, GSVA analysis was employed to calculate the enrichment scores of HALLMARK and KEGG pathways. Interestingly, some significantly enriched HALLMARK pathways that may be directly related to the dilation of the heart of DCM patients were detected.
For example, heme metabolism was highly enriched in subtype-1, and myogenesis was enriched in subtype-2 ( Figure 6A). Besides, several cancer related KEGG pathways were dramatically enriched in subtype-1, such as renal cell carcinoma, pancreatic cancer, prostate cancer and melanoma ( Figure 6B). In contrast, various metabolic

| Co-expression network analysis identified autophagy expression pattern related gene modules
Next, a comprehensive gene landscape correlated to each autophagy expression patterns was constructed, and gene coexpression modules related to distinct autophagy regulations were identified by weighted gene co-expression network analysis (WGCNA) (Figure 7A,B). Seven gene modules were determined and different expression pattern matched their related genes ( Figure 7C, Table S15). Among them, the most relevant to the autophagy subtype was the blue module ( Figure 7D), which was closely correlated to the subtype-2 (cor = 0.99, p < 1e-200).
In order to further figure out the regulatory roles of autophagy genes in the blue module, the interaction network including DCMrelated autophagy genes and the hub genes (kME >0.95) in this F I G U R E 5 Diversity of immune microenvironment characteristics between distinct autophagy-mediated regulation patterns. (A) The abundance differences of infiltrating immunocytes between 2 autophagy regulation patterns. (B) The activity differences of immune reaction gene-sets between 2 autophagy regulation patterns. (C) The expression differences of HLA genes between 2 autophagy regulation patterns ns ns **** * **** ns ns **** ns *** ns **** ns * ns * * ns ns ****

| DISCUSS ION
Autophagy is a bridge connecting innate immunity and adaptive immunity. 21  Recent years, although gene sequencing technology based on the gene mutation analysis has been used for the diagnosis of DCM, 5,7,48 this relies heavily on the detailed information of family genetic history. Besides, the accuracy of this approach needs to be improved.
With the outbreak of COVID-19, the nucleic acid detection method has become the most popular genetic diagnosis technology. In this study, 13 DCM-related autophagy genes were screened through a series of bioinformatics methods and a risk diagnosis model was constructed, which can distinguish healthy and DCM samples well ZNF358 sex hormones to their receptors directly changes the functions of immune cells and platelets, thereby affecting the type of cardiac inflammation, remodelling and thrombosis in DCM. Interestingly, the androgen response pathway was highly enriched in autophagy subtype-1, while the oestrogen response pathway was highly enriched in subtype-2. This indicated that these two distinct autophagy-mediated expression patterns respond to different sex hormones, which may lead to gender differences in the risk of DCM. In addition, some cancer-related pathways were highly enriched in subtype-1, and the KEGG enrichment analysis of the differential expressed genes between the two subtypes also showed that the most significantly enriched pathways were related to diseases. This suggested that the regulation of autophagy may be involved in the complications of DCM and other diseases.
Overall, our findings revealed how autophagy affects the immune microenvironment of DCM, and provided new insights to understand the pathogenesis of DCM. This study was the first to systematically reveal the potential connections between autophagy and immune microenvironment in DCM. These findings provide clues for further studying the mechanism of autophagy in DCM.

CO N FLI C T O F I NTE R E S T
The authors declare that they have no competing interests.

DATA AVA I L A B I L I T Y S TAT E M E N T
The data that support the findings of this study are available in