Shangjie Zou1,2
Bee Luan Khoo1,2,3*- 1Department of Biomedical Engineering, City University of Hong Kong, Hong Kong, Hong Kong SAR, China
- 2Hong Kong Center for Cerebro-Cardiovascular Health Engineering (COCHE), Hong Kong, Hong Kong SAR, China
- 3Department of Precision Diagnostic and Therapeutic Technology, City University of Hong Kong-Shenzhen Futian Research Institute, Shenzhen, China
Background: A central issue hindering the development of effective anti-fibrosis drugs for heart failure is the unclear interrelationship between fibrosis and the immune cells. This study aims at providing precise subtyping of heart failure based on immune cell fractions, elaborating their differences in fibrotic mechanisms, and proposing a biomarker panel for evaluating intrinsic features of patients’ physiological statuses through subtype classification, thereby promoting the precision medicine for cardiac fibrosis.
Methods: We inferred immune cell type abundance of the ventricular samples by a computational method (CIBERSORTx) based on ventricular tissue samples from 103 patients with heart failure, and applied K-means clustering to divide patients into two subtypes based on their immune cell type abundance. We also designed a novel analytic strategy: Large-Scale Functional Score and Association Analysis (LAFSAA), to study fibrotic mechanisms in the two subtypes.
Results: Two subtypes of immune cell fractions: pro-inflammatory and pro-remodeling subtypes, were identified. LAFSAA identified 11 subtype-specific pro-fibrotic functional gene sets as the basis for personalised targeted treatments. Based on feature selection, a 30-gene biomarker panel (ImmunCard30) established for diagnosing patient subtypes achieved high classification performance, with the area under the receiver operator characteristic curve corresponding to 0.954 and 0.803 for the discovery and validation sets, respectively.
Conclusion: Patients with the two subtypes of cardiac immune cell fractions were likely having different fibrotic mechanisms. Patients’ subtypes can be predicted based on the ImmunCard30 biomarker panel. We envision that our unique stratification strategy revealed in this study will unravel advance diagnostic techniques for personalised anti-fibrotic therapy.
Introduction
Heart failure results from the progression of cardiovascular diseases, such as coronary artery diseases and cardiomyopathy (1). There are two types of heart failure: ischemic heart failure and non-ischemic heart failure. Ischemic heart failure is typically caused by blood vessel diseases that result in myocardial infarction and cardiomyocyte death (1), whereas non-ischemic heart failure is caused by genetic mutations, infections, or exposure to stimuli that result in myocardium dysfunction (2). Early-stage heart failure is typically accompanied by loss of functional cardiomyocytes, hypertension, and increased ventricular wall stress, all of which contribute to inflammation and pro-fibrotic signaling (3). As a result, fibrosis is a common pathological process in the progression of heart failure caused by both etiologies (2, 3), and it is commonly regarded as a biomarker for poor prognosis (4–6) as well as the hallmark of end-stage heart failure (7).
Although fibrosis was initiated as a means of healing, it would eventually lead to decreased ventricular compliance, abnormal cardiac rhythm, and cardiomyocyte death, worsening the symptoms of heart failure (8). Conventional heart failure therapies, which focuses on symptom relief rather than maladaptive remodeling processes such as cardiac fibrosis (9, 10), would lead to varied outcomes among patients (11). Current clinical research on anti-fibrosis treatments has focused on directly reversing the activated cardiac fibroblast phenotype (11, 12). However, most of these treatments remain unsatisfactory.
Previous research has shown that, while activation of immune responses and immune cells is required for cardiac repair, excessive activation is a major cause of fibrosis (3, 13). Nonetheless, the lack of understanding of the interplay between cardiac fibroblasts and the immune cell fractions prevents therapeutic agents from targeting the pro-fibrotic regulatory networks underlying fibroblast activation (12, 14–18). For example, endothelin inhibitors, anti-TGFβ, and anti-inflammatory drugs target fibrotic pathways but had a high degree of inconsistency in treatment efficacy within patients (19, 20), and were even reported to increase the risk of adverse immune responses (21). These unsatisfactory results can partially due to our limited understanding of the complex relationships between cardiac fibrosis and immune cells. For example, the overlapping impact of M1 and M2 macrophages on promoting fibrosis through different approaches is attributed to the dual-pronged roles of immune activation, despite having opposite regulatory functions on inflammation (13). Neutrophils, monocytes, T cells, and B cells were thought to regulate fibrosis, but crosstalk and fibrosis regulation (promoting or suppressing) are related to pathological statuses and immune cell populations (13, 22). CD4+ T cells mainly promote fibrosis through secreting cytokines like IFN-γ, IL-4, IL-13, and IL-17 (23). B cells can secret immunoglobulins to directly induce fibrosis (24). Monocyte promotes inflammation and immune cell infiltration (13). M2 Macrophage reduce inflammation, but promotes fibrosis through the process of tissue healing and cardiac remodeling (13). These evidences in former reports indicates that immune cell fraction is an important causal factor for fibrotic mechanisms in heart failure.
As a result, finding a therapy targeting immune response that can achieve consistent results in all patients would be difficult without considering their local conditions of immunity on a more systematic scale.
Understanding the interrelationship between cardiac fibroblasts and their immune cell fractions holds promise for breakthroughs in effective anti-fibrosis therapy (11). Furthermore, heart failure and cardiac fibrosis are also present with complex symptoms and multiple etiologies, which brings a tremendous amount of challenge to achieving individualised treatment strategies (25, 26). There is an urgent need for a better understanding of the heterogeneity in immune cell fractions and its relationship to fibrosis as an important prerequisite for effective anti-fibrosis strategies.
Here, we report for the first time the establishment of a matrix of cardiac immune cell fractions for heart failure patients and validated the subtyping analysis based on that. Based on cardiac immune cell fractions, two distinct patient subtypes were identified: [1] pro-inflammatory (with elevated levels of biomarkers related to adaptive immune responses) and [2] pro-remodeling (with high M2 macrophage). We demonstrated that epithelial-mesenchymal transition (EMT) was highly correlated with fibrosis in both subtypes. A novel analytical method based on Large-Scale Functional Score and Association Analysis (LAFSAA) identified functional gene sets related to subtype-specific mechanisms of EMT and fibrosis. Fibrosis of the pro-inflammatory subtype was specifically correlated with chronic inflammatory responses, JAK activities, and hemidesmosome. In contrast, fibrosis of the pro-remodeling subtype was specifically correlated with IL-6 receptor binding, glycoside metabolic process, leukocyte aggregation, and blood vessel development. We then present a unique biomarker panel, termed ImmunCard30, based on 30 genes encoding secretory proteins for patient stratification. The signature score of ImmunCard30 achieved high subtype classification efficiency, with the area under the receiver operator characteristic curve (AUC) equal to 0.958 in the discovery set (n = 103) and 0.79 in the validating set (n = 177). ImmunCard30 could potentially serve as the biomarker panel for diagnosing patients’ subtypes and guide therapeutic treatments targeting the subtype-specific pro-fibrotic mechanisms. Overall, this study provides new insights into the heterogeneous cardiac immune cell fractions and its interactions with fibrosis, thereby providing an important basis for anti-fibrotic drug screening and disease management strategies. The novel patient stratifying strategy and the ImmunCard30 biomarker panel provide the basis for an important step towards precise anti-fibrotic therapies that surpass the current treatments.
Materials and methods
Data installation and preparation
The transcriptomic datasets were obtained from the Gene Expression Omnibus (GEO) database. Five datasets of ventricular tissue samples were involved, among which three of them were used to form the discovery set: GSE116250 (27), GSE135055 (28), and GSE46224 (29). Independent dataset GSE57338 (30) was used as the validating set. In total, 103 samples were included in the discovery set, and 177 samples were included in the validating sets. Two datasets from lung tissues were involved in validating signature scores of fibrosis: GSE124685 (31) and GSE184316. Raw data of RNA-sequencing was transformed into reading counts using Genome Reference Consortium Human Build 38 as the reference genome. The read counts of RNA-sequencing were further normalised to Reads Per Kilobase of transcript per Million reads mapped (RPKM) by R package edgeR (R 4.1.0) (32).
Apart from the transcriptomic datasets, this study also involved information from public biomedical databases. Annotated gene sets of EMT and angiogenesis were retrieved from hallmark gene sets of MSigDB collections (33, 34). The list of secretory proteins was obtained from The Human Protein Atlas (https://www.proteinatlas.org/) (35).
Quantification and clustering of the immune cell fractions
CIBERSORTx calculated the matrix of samples’ immune cell fractions based on the 22-leukocyte signature matrix (LM22) (36, 37). The calculated matrix of immune cell fractions serves as the input of clustering to identify subtypes of the immune cell fractions. The average silhouette width was applied to evaluate optimal numbers of clusters (38), while k-means clustering was conducted based on R package factoextra (R 4.1.0) (39). Visualisation of k-means clustering was achieved by dimensional reduction based on t-distributed stochastic neighbor embedding (tSNE) implemented by R package RtSNE (R 4.1.0) (40). Recursive feature elimination for ranking immune cell types on their importance in determining subtype affiliations was achieved by R package caret (41).
Transcriptomic analysis and statistical analysis
Differential gene expression analysis was conducted based on R package DESeq2 (R 4.1.0) based on read counts of RNA-sequencing normalised by datasets (to remove potential batch effects) (42). Functional enrichment analysis was conducted by the online tool Metascape (https://metascape.org/) to study enriched gene sets in Gene Ontology Biological Process, Kyoto Encyclopedia of Genes and Genomes, and Hallmark gene sets (43). Correlations between variables were evaluated by Pearson correlation with R scripts (4.1.0). Signature scores were calculated based on R package singscore (R 4.1.0) (44). Scores of fibrosis, scores of EMT, and scores of functional gene sets were calculated by singscore using a uni-directional mode.
In contrast, the score of the predictive gene panel for subtype classification was calculated based on the bi-directional singscore method (“bi-directional” means considering both positively- and negatively-correlated genes in the scoring). Feature extraction based on least absolute shrinkage and selection operator (LASSO) regression was performed by R package glmnet (R 4.1.0), using the variables corresponding to the highest ln (λ) value that resulted in a misclassification error within 1 standard error of the minimum cross-validated error as the optimal variable combination (45). ROC curves and AUC were plotted and calculated by R package pROC (R 4.1.0) (46).
Implementation of LAFSAA
By applying singscore (44) method, we calculated signature scores of the 10402 ontology gene sets and 50 HALLMARK gene sets retrieved from MSigDB (33) based on the expression matrix of the 103 samples from the discovery set. The matrix containing samples on the rows, functional gene sets on the columns, and signature scores as values was named the quantitative functional matrix. Spearman correlation analysis was applied to screen the gene sets in the quantitative functional matrix. Those gene sets that were significantly and specifically correlated with the target gene set (TGS, which was the HALLMARK gene set of EMT in this study) in samples from one subtype would be identified as pro-inflammatory or pro-remodeling subtype-specific gene sets. The identified gene sets would be further screened based on their correlations with EMT in their corresponding subtypes in the validating set GSE57338.
In the meantime, to make the screening more stringent, we applied concept signature enrichment analysis (CSEA) (47) based on the genes that were differentially expressed between the high and low fibrosis groups in each subtype. Only gene sets with p-value< 0.05 in CSEA analysis of their corresponding subtypes were regarded as subtype-specific EMT-related gene sets.
Calculation of intra-cluster correlation scores
Intra-cluster correlation scores were used to select representative gene sets for clusters identified based on the correlation matrix. It denotes the average coefficients of a gene set’s correlations with other gene sets in its corresponding cluster. This score was calculated by the equation below:
Where ICS(GSi) is the intra-cluster correlation score of the ith gene set. The n denotes the number of gene sets in a cluster. SpearmanCor(GSi, GSj) denotes the spearman correlation of the ith and the jth gene set.
Results
Clustering of immune cell fractions for subtyping
Samples from patients with end-stage heart failure (n = 103) and healthy samples (n = 31) from three datasets were retrieved (Table 1). Since infiltrated immune cells are the major components and central regulators of cardiac immune responses, we applied CIBERSORTx to calculate immune cell fractions for all the involved patient samples (36) (Supplementary Table 1), and conducted subtyping based on the matrix of immune cell fractions.
Table 1 Clinical characteristics of samples from discovery set.
We identified two potential clusters based on the matrix of immune cell fractions, which was the average silhouette width method (Supplementary Figure 1). We demonstrated that k-means clustering classified the tissue samples (n = 103) into two subtypes, termed Subtype 1 and 2 (Figure 1A). Samples were assigned to Subtype 1 (n = 34) and Subtype 2 (n = 69), respectively (Figure 1A). Checking on samples’ etiologies and source datasets indicated that clustering results were not influenced by dataset attribution or underlying diagnosis (Figure 1A; Supplementary Table 2).
Figure 1 Clustering of heart failure samples for subtyping stratification. (A) Left panel: clustering of the two subtypes was visualized by t-distributed stochastic neighbor embedding (t-SNE). Red dots represented subtype 1 samples, while turquoise dots represented subtype 2 samples. Right panel: distribution of heart failure samples from different datasets (upper panel) and with different diagnoses (lower panel) in the two subtypes. Datasets were displayed with their respective Gene Expression Omnibus Series (GSE) identities. The samples’ diagnosis was categorized into ischemic cardiomyopathy (ICM) and non-ischemic cardiomyopathy (NICM). (B) Boxplot showing the immune cell fractions of the 22 cell types in the two subtypes. Red and cyan boxes, respectively, represented subtypes 1 and 2. The significance of the difference in immune cell fractions between the two subtypes calculated by the two-tailed Student’s t-test was shown on the top of each cell type’s column. (C) Heatmap showing the matrix of significantly different immune cell fractions between the two subtypes. Rows represented the immune cell types annotated on the right side of the heatmap. Each column represented one sample, while column annotations on the top of the heatmap described the corresponding samples’ subtype belongings and fibrotic status. (D) Boxplot showing the signature score of fibrosis in the two subtypes. Red and cyan boxes, respectively, represented subtypes 1 and 2. Each black point represented one sample in its corresponding subtype. The significance of the difference in score between the subtype 1 (n = 34) and the subtype 2 (n = 69) calculated by the two-tailed Student’s t-test (t = 2.187, df = 101) was shown at the top of the plot. For each boxplot: the center line represents the median, the box represents the interquartile range, and the whisker displays minimum to maximum. (P-value< 0.05: *< 0.01: **< 0.001: ***< 0.0001: ****). (E) Scatter plots showing log2(fold-change) of genes in the two subtypes compared with healthy controls. Each dot represented one gene. In the upper panel, genes with log2(fold-change) higher than 1 or lower than -1 in one of the subtypes were highlighted in red or blue, respectively. In the lower panel, genes with log2(fold-change) in subtype 1 that were two times higher than in subtype 2 were highlighted in red. (F) Functional enrichment analysis of the genes highlighted in red in (E) lower panel. The bar plots were -log10 (P-value) of enrichments of genes in signatures of Gene Ontology (GO) Biological Process and Kyoto Encyclopedia of Genes and Genomes (KEGG) of the Molecular Signatures Database.
Based on boxplots showing comparisons of immune cell proportions, type 2 subtypes were higher in neutrophils and M2 macrophages, whereas type 1 was higher in B cells, T cells, and monocytes (Figure 1B). Figure 1C shows a heat map showing significantly different immune cell levels between the two subtypes.
To quantify the fibrotic status of the samples, we need to select marker genes that can represent fibrosis. POSTN, COL1A1, TIMP1, MMP2, and ACTA2. POSTN and COL1A1 are genes that express the extracellular matrix proteins periostin and type 1 collagen, both of which are significantly up-regulated in heart failure with cardiac fibrosis (48). TIMP1 and MMP2 are enzymes that aid in the progression of cardiac fibrosis (49). ACTA2 is a fibroblast activation biomarker (50). In previous studies, these direct fibrosis marker genes were used to assess pro-fibrotic activity (51). The same five-marker combination has been used to indicate fibrotic levels in patients (50). Fibrosis signature scores based on these five biomarkers were validated on multiple organ fibrosis datasets corresponding to heart and lung tissue, showing significant differences between normal and fibrotic samples (Supplementary Figure 2). We applied signature scoring to the 103 tissue samples and divided the samples into high and low fibrosis groups based on the median of fibrotic signature scores. Overall, Subtype 1 (n = 34) had a slightly higher fibrosis score than Subtype 2 (n = 69) (Figure 1D), but both high and low fibrosis samples were seen in both subtypes (Figure 1B).
To further evaluate the differences in fibrosis-related factors between the two subtypes at the molecular level, we explored the differentially expressed genes (DEG) separately in the two subtypes relative to healthy controls. DEGs in both subtypes had similar fold changes to healthy controls (Figure 1E). However, 286 genes showed enhanced heart failure-related up-regulation in Subtype 1 (2 times higher in log2 fold change) than Subtype 2 (Supplementary Table 3; Figure 1E).
Among the genes with the highest levels of upregulation in subtype 1: The functional enrichment analysis ranked the Gene Ontology gene set “adaptive immune response” first, with multiple other gene sets related to adaptive immune responses (for example, regulation of lymphocyte activation, T cell selection, T cell receptor signaling, and regulation of B cell proliferation) identified as significantly enriched as well (Figure 1F). Other inflammatory signaling pathways involved in adaptive immune responses, including those related to IL-2, JAK-STAT signaling, and IL-17 signaling, were also more activated in the subtype 1 (Figure 1F).
As a further demonstration, we selected a list of genes as typical biomarkers of adaptive immune response, which included multiple surface antigens of T cells (such as CD8A, CD8B, CD3G), NK cells (such as CXCR3, CD7), or B cells (such as CXCR5, CD79A), and genes involved in signaling of chemokines and cytokines that regulates adaptive immune responses. According to the differential gene expression analysis, these genes showed higher up-regulation compared to healthy controls in the subtype 1 than in the subtype 2 (Supplementary Figure 3A).
In the meantime, the subtype 1 was higher in proportion of CD4+ T cells, regulatory T cells, memory B cells, and plasma cells (Figure 1B). The high percentage of T and B cells indicated an active adaptive immune response. Furthermore, to directly compare the differences in activation of T cells and B cells between the two subtypes, we calculated signature scores of multiple functional gene sets related to B cell activation, T cell proliferation, and T cell extravasation (all of which are crucial processes in adaptive immune responses), and demonstrated that subtype 1 was higher in all of these gene sets (Supplementary Figure 3B). Meanwhile, so as to evaluate the activity of B cells, we have also calculated the expressional score of immunoglobulins based on expression level of 188 genes that express immunoglobulins (calculated by “singscore” method). The obtained expressional score of immunoglobulins was significantly higher in the subtype 1, indicating a higher activity of antibodies in the subtype 1 (Supplementary Figure 3B).
In conclusion, biomarker and functional gene set analyses revealed that the two immune cell subtypes were associated with different levels of biomarkers related to adaptive immune responses. Further research on patient samples using comprehensive clinical recordings or rationally designed disease models is expected to provide more evidence to validate this correlation in the future.
Distinct immune cell types related to subtype stratification and fibrosis
We scored and ranked 22 leukocyte cell types based on recursive feature elimination according to their importance in determining subtype affiliations (Figure 2A). To validate the ranking, we applied CIBERSORTx to an independent dataset with a large sample size (n = 177): GSE57338 (30) and assigned samples to two subtypes by calculating their distance from the centroids of the two subtypes in the discovery set (Supplementary Figure 4).
Figure 2 Immune cell types related to subtype stratification and fibrosis. (A) Rank of immune cell types based on their scores calculated by recursive feature elimination. (B) Pearson correlations between immune cell types. The colored square in each rectangular cell displayed the correlation coefficient between immune cell types on its corresponding row and column. Red squares indicated negative coefficients, while blue squares indicated positive coefficients. Darker shades and larger square sizes indicated a higher absolute value of coefficients. Only correlations with P-value less than 0.05 were displayed with correlation coefficients. (C) Scatter plots with Pearson correlation coefficients between fractions of each immune cell type and the signature score of fibrosis. Blue light in each subgraph represented the linear regression line. Each black dot represented a heart failure sample. The light gray area surrounding each blue regression line represented a level 0.95 confidence interval.
According to the result of recursive feature elimination, the proportion of M2 macrophages, which was significantly higher in the subtype 2, was robustly the most important feature that discriminated the two subtypes (Figure 2A). The proportion of resting CD4+ T cells was identified as the most important feature for the subtype 1 (Figure 2A). With the distinct phenotype of elevation in resting CD4+ T cell proportions, leading to high immune response activation, Subtype 1 was subsequently referred to as the pro-inflammatory subtype. The high level of M2 macrophages, which plays a pivotal role in anti-inflammation and pro-remodeling in fibrotic tissues (52), is the distinct phenotype of Subtype 2, termed the pro-remodeling subtype.
Figure 2B depicts a broad overview of immune cell type correlations. CD4+ T cell levels were found to be positively correlated with monocytes, M0 macrophages, and memory B cells, indicating that these cell types have close regulatory relationships in the pro-inflammatory subtype. M2 macrophage levels, on the other hand, were significantly inversely correlated with the proportions of monocytes, T cells, and B cells (Figure 2B), confirming the role of macrophage M2 in inflammation suppression (53) and its potentially key role in determining the immunity of the pro-remodeling subtype, which was distinct from the pro-inflammatory subtype.
We calculated the Pearson correlation coefficient between the fractions of various immune cells and fibrosis to assess the potential contribution of immune cells to fibrosis (Figure 2C). The immune cells chosen were either important for classifying the two subtypes (for example, M2 macrophages, CD4+ resting T cells, memory B cells, and M0 macrophages) or had high proportions in the immune cell population (such as NK cells, mast cells). The proportion of monocytes was found to be positively correlated with fibrosis. The proportion of neutrophils, on the other hand, was negatively correlated. Despite being the most distinguishing features of the two subtypes, the proportions of M2 macrophages and CD4+ T cells were not significantly associated with fibrosis in our discovery set (Figure 2C).
EMT was potentially a major contributor to fibrosis in both subtypes
In the previous analyses, we demonstrated that the major features of the two subtypes: the proportion of M2 macrophages and CD4+ T cells, were not significantly related to fibrosis (Figure 2C). As a result, subtype affiliation may not be a direct predictor of fibrotic levels. To better understand the fibrotic mechanisms in these two subtypes, we evaluated if there was a universal hub pro-fibrotic process induced in both subtypes.
By dividing the tissue samples (n = 103) using the median fibrosis score as a cutoff, we obtained high and low fibrosis groups in both the pro-inflammatory subtype (High: n = 20, 19.42%; Low: n = 14, 13.59%) and the pro-remodeling subtype (High: n = 31, 30.10%; Low: n = 38; 36.89%). The distribution of high and low fibrosis samples in the two subtypes can be seen in the column annotation of the heatmap Figure 1C. Differential gene expression analysis between the high and low fibrosis groups identified 166 commonly up-regulated genes in the high-fibrosis group of both subtypes (Supplementary Table 4; Figure 3A).
Figure 3 Identification of epithelial-mesenchymal transition (EMT) as the potentially pivotal pro-fibrotic process. (A) Expression of 166 common genes that were up-regulated in high fibrosis samples in both subtypes. The expression level was displayed by z-score. Each yellow line represented the line plot of a gene’s expression level across the 103 heart failure samples and 31 healthy samples. The black line represented the average of z-scores of the 166 genes among the samples. Vertical dash lines separated the samples from different groups, with their belongings shown below x-axis. (B) Bubble plot for functional analysis of the 166 up-regulated genes in high fibrosis samples. The deeper red color of the bubbles represented a higher -log10(P-value), while the larger circles indicated the higher number of genes identified in the corresponding gene set. (C) Pearson correlation coefficients between signature scores of fibrosis, EMT, and angiogenesis. Numbers and asterisks in rectangular cells represented correlation coefficients and P-values (ns, nonsignificance; Pvalue >0.05: *> 0.01: **> 0.001: ***> 0.0001: ****). (D) Pearson correlation coefficients (left panel) and -log10(P-value) of the correlations (right panel) between expression levels of the 46 key EMT-related transcription factors and score of fibrosis. Levels of correlation coefficients and -log10(P-value) were displayed with colors according to the legends on the right side of their respective panel. Only those genes with -log10 (P-value) > 1.30103 (P-value< 0.05) were displayed with color in the right panel.
Functional enrichment analysis identified that a large number of commonly up-regulated genes in the high fibrosis group were associated with EMT (Figure 3B). Signature scores calculated based on the HALLMARK gene set of EMT demonstrated that EMT was highly up-regulated in heart failure patients, and the high fibrosis groups had significantly higher scores than low fibrosis groups (Supplementary Figure 5). There are normally two types of EMT in adult ventricular tissues: epicardial EMT and endothelial-mesenchymal transition (EndMT) (54). EndMT is commonly regarded as a subset of EMT due to shared pathways and transcriptional factors (55). The signature scores of EMT and angiogenic hallmark signatures were highly correlated with each other (coefficient = 0.73), as was the signature score of fibrosis (coefficient equals 0.86 for correlation between EMT and fibrosis, and equals 0.62 for correlation between angiogenesis and fibrosis) (Figure 3C), indicating potential relationships between vascular endothelial cell activities and EMT in the discovery set.
We studied the correlations of cytokeratin (expressional score calculated by “singscore” method based on 55 genes expressing cytokeratin) and PECAM1 (CD31) with EMT and fibrosis to determine which type of EMT was the major contributor to fibrosis. Previously, cytokeratin was thought to be an epithelial cell feature, while PECAM1 was thought to be an endothelial cell feature (56, 57). The correlation coefficient between PECAM1 and EMT was 0.42, with a high significance (P 0.001), while the correlation between cytokeratin and EMT was much lower (coefficient = 0.24, P 0.05) (Figure 3C). These results showed that EMT in the discovery set was closely related to endothelial cells.
We discovered that endothelial biomarkers were highly correlated with EMT and fibrosis in patients with end-stage heart failure. Previous research has suggested that EndMT may contribute to cardiac fibrosis (58). As a result, more clinical cohort studies may be needed to confirm whether the identified EMT in end-stage heart failure patients, which was highly correlated with fibrosis, was related to EndMT.
At the level of transcriptional regulation, we screened 46 essential transcriptional factors associated with EMT based on previously reported regulatory landscapes (Supplementary Table 5) (59), of which 29 were significantly associated with fibrosis scores in the discovery set (all significant correlation coefficients were positive) (Figure 3D). These significant correlations served as strong evidences to prove the associations between EMT-related transcription factors and fibrosis, validating the pivotal role of EMT as a central process in fibrosis.
Mechanistic study of EMT in the two subtypes based on LAFSAA
The stratification of heart failure patients into pro-inflammatory and pro-remodeling subtypes would unravel unique biomarkers related to their fibrotic mechanisms. Here we designed a novel approach based on LAFSAA to identify the functional gene sets that may be associated with distinct EMT mechanisms in the two subtypes (Supplementary Figure 6). We identified 17 gene sets specifically associated with EMT in the pro-inflammatory subtype and 87 gene sets specifically associated with EMT in the pro-remodeling subtype.
To study the potential functional relationships between these 104 gene sets, we calculated their correlation matrix based on Spearman correlation coefficients in the discovery set. According to the correlation matrix, hierarchical clustering divided the 104 gene sets into 8 clusters with high intra-cluster correlations (Figure 4A). Gene sets assigned to the same cluster were expected to be functionally relevant. Gene sets with intra-cluster correlation scores (Calculated by Eq. 1) ranked in the top 30% of their corresponding clusters were especially noted to interpret functions of their corresponding clusters (Figure 4A).
Figure 4 LAFSAA-identified subtype-specific pro-EMT and pro-fibrotic gene sets. (A) Heatmap on the top left shows the correlation matrix of LAFSAA-identified subtype-specific pro-EMT gene sets. Column annotation of the heatmap indicates the subtype-specificity of these pro-EMT gene sets. The gene sets were grouped into 8 clusters based on hierarchical clustering. Intra-cluster correlations were displayed in the form of a correlation network, with edges denoting the correlation coefficients between each pair of functional gene sets. Gene sets that were represent the functions of their corresponding clusters were marked in magenta, with their intra-cluster correlation scores displayed in their corresponding bar charts. (B) Subtype-specific pro-fibrotic gene sets and their Spearman correlation coefficients with fibrosis in the two subtypes in the discovery set and the validating set GSE57338. The gene sets were marked with their belonged LAFSAA-identified clusters.
Among the LAFSAA-identified functional clusters, clusters 1 (related to glycoside metabolic process), 3 (related to the immune response to environmental stresses), 4 (related to rhythmic process and neural development), 5 (hyaluronan biosynthesis and metabolism), and 7 (related to cell adhesion, cytoskeletal organisation, and vasculature development) were specifically correlated with EMT in the pro-remodeling subtype. In contrast, cluster 2 (related to translational regulation, neural system development, and cell adhesion) and 6 (related to JAK activities, inflammation and adaptive immune response) were specifically correlated with EMT in the pro-inflammatory subtype (Figure 4A).
Among the subtype-specific EMT-correlated gene sets from these 8 clusters, 11 of them (from clusters 1, 2, 3, 6, and 7 of LAFSAA-identified functional clusters) were also subtype-specific fibrosis-correlated gene sets, as shown in Figure 4B. Meanwhile, these 11 functional gene sets were mainly having similar signature scores between the two subtypes (Supplementary Figure 7). These results indicated that the 11 functional gene sets were likely represented the heterogeneous fibrotic mechanisms in different subtypes of immune cell fractions.
Among these 11 gene sets, 7 were specifically correlated with fibrosis in the pro-inflammatory subtype (Figure 4B). These 7 gene sets were mainly related to immune responses, including “GOBP: regulation of chronic inflammatory response” and “GOBP: activation of Janus kinase activity”. This result demonstrated that in the pro-inflammatory subtype, fibrosis was potentially promoted by EMT induced by inflammatory responses. Meanwhile, hemidesmosome, regulation of protein localization to synapse, and limb bud formation (related to cell morphogenesis) were also specifically correlated with EMT and fibrosis in the pro-inflammatory subtype, which might reflect the reorganization of extracellular matrix and morphological changes of cardiac cells.
In contrast, gene sets specifically correlated with fibrosis and EMT in the pro-remodeling subtype were related to IL-6 receptor binding, glycoside metabolism, leukocyte aggregation, and varicosity. IL-6 was identified as a key factor regulating vascular remodeling in former studies (60). In this study, IL-6 receptor binding and varicosity were assigned to the same cluster in the similarity matrix of LAFSAA (Figure 4A), which indicated that IL-6 and vascular remodeling were at least statistically relevant in the samples of the discovery set, and potentially being functionally relevant in the process of EMT and progression of fibrosis. In former studies, M2 macrophage was believed to strongly promote vascular development (61). By combining the analytical results in this study and the physiological relationships between M2 macrophage and vascular development revealed in former studies, the high correlation between vascular remodeling and fibrosis in the pro-remodeling subtype is potentially due to the high proportion of M2 macrophage in immune cell population.
We proposed functional gene sets with subtype-specific correlations with EMT and fibrosis via LAFSAA. These findings suggested that the two subtypes had distinct fibrotic mechanisms. More research using larger clinical datasets or disease models representing the two subtypes of immune cell fractions are needed for further validating these correlations. More importantly, the different fibrotic mechanisms of the two subtypes heralds that patient stratification based on subtypes of immune cell fractions could potentially benefit anti-fibrotic therapies. Inhibiting the activities of subtype-specific fibrosis-related functional gene sets based on patients’ subtype affiliations could be a promising approach to improving therapeutic outcomes compared to standard therapeutic strategies, which should be investigated further with in-vivo and in-vitro tests.
Patient stratification using immune subtyping for personalised therapy
The key benefit of precise stratification of patients is the ability to target fibrotic mechanisms associated with the immune cell fractions for effective treatment. The design of a strategy to predict subtypes of cardiac immune cell fractions in patients and assist in the development of treatment strategies is pivotal. However, current subtyping based on RNA profiling data from ventricular tissue samples cannot be directly applied in clinical settings. Therefore, we aimed to identify a panel of biomarkers correlated with subtypes of cardiac immune cell fractions that could be broadly applied for practical applications in personalized medical diagnosis.
We selected genes encoding secretory proteins for subtyping as potential biomarkers for assay-based screening. The Student’s t-test identified genes encoding secretory proteins that were expressed significantly different between the two subtypes (Supplementary Figure 8). The differentially expressed genes were screened with the LASSO regression (Supplementary Figure 9). We identified a unique panel of 32 genes with high performance in predicting patients’ subtypes of cardiac immune cell fractions using LASSO screening (Supplementary Figure 10A). To improve the robustness of the biomarker panel, we applied receiver operating characteristic (ROC) analysis to evaluate the classifying efficiency of genes that were not highly significantly different between the two subtypes according to the heatmap shown in Supplementary Figure 10A. We demonstrated that excluding HLA-DPA1 or KRT9 from the panel led to an increase in classifying efficiency in the validating set GSE57338, but not to a significant decrease in AUC in the discovery set, by successively excluding the selected genes from the panel and recording the changes in classifying efficiency quantified by AUC (Supplementary Figure 10B).
After excluding HLA-DPA1 and KRT9, the AUC of the biomarker panel in validating set GSE57338 subtype classification increased from 0.79 to 0.802, while the AUC in the discovery set decreased from 0.958 to 0.954. (Supplementary Figure 10B). As a result of this finding, HLA-DPA1 and KRT9 were removed from the biomarker panel. Therefore, the optimal combination of biomarker panel consisted of 30 genes, termed ImmunCard30. In the ImmunCard30 panel, 17 genes were expressed significantly higher in the pro-inflammatory subtype (hereafter termed as 17G-subtype1), and 13 genes were expressed significantly higher in the pro-remodeling subtype (hereafter termed as 13G-subtype2) (Figure 5A)
Figure 5 Performance evaluation for ImmunCard30 in predicting patients’ subtypes. (A) The average expression level of the 30 genes for subtype classification. The average z-score of each gene in the corresponding subtype was used to represent the expression level. (B) Signature scores of the two subtypes in the discovery set. (C) Signature scores of the two subtypes in the validating set GSE57338 were compared. (D) Receiver operating characteristic curves of the 30-gene biomarker panel was used to classify the two subtypes from the discovery set. (E) Receiver operating characteristic curves of the 30-gene biomarker panel was used to classify the two subtypes from the validating set GSE57338. A two-tailed Student’s t-test was used to determine significance. AUC, area under the curve.
To validate the performance of ImmunCard30, we applied the bi-directional singscore method to calculate signature scores for the clinical samples, where samples with a higher expression level of the 17G-subtype1 would be given higher scores, and a higher expression level of the 13G-subtype2 would lead to lower scores. We demonstrated that the samples from the pro-inflammatory subtype had significantly higher signature scores of the ImmunCard30 predictive panel than those from the pro-remodeling subtype, as confirmed in both the discovery and the validating set GSE57338 (Figures 5B, C). Signature scores calculated by the ImmunCard30 predictive panel achieved a high AUC of 0.954 in the discovery set (n = 103) and a relatively good predicting result in the validating set’s 177 samples (AUC = 0.803) (Figures 5D, E).
The functional correlations between ImmunCard30 and fibrotic mechanisms could possibly be an approach to better understand the different fibrotic mechanisms in the two subtypes. MSigDB database contains over 10000 functional gene sets, with clear recordings of their functional annotations and gene lists. We screened the functional gene sets in the MSigDB database and identified those fibrosis-correlated functional gene sets (with signature scores significantly correlated with fibrosis, with a p-value of Spearman’s correlation< 0.05) that contained genes within the ImmunCard30. Functional annotations of these fibrosis-correlated functional gene sets from MSigDB were used to interpret the functions of genes from the ImmunCard30. According to the results, we found that genes within the ImmunCard30 panel were related to multiple pro-fibrosis processes or molecular components, such as adaptive immune response, cell migration and adhesion, neural system development, and vasculature development (Figure 6). In general, biomarkers in the 17G-subtype1 (ImmunCard30 genes that expressed significantly higher in the pro-inflammatory subtype) were more related to adaptive immune response, inflammation, and neural system development. In contrast, biomarkers in the 13G-subtype2 (ImmunCard30 genes that expressed significantly higher in the pro-remodeling subtype) were specifically related to vasculature development (Figure 6). Biomarkers from both subtypes were related to immune responses, cell migration, and cell adhesion (Figure 6).
Figure 6 Functional annotations of the 30 biomarkers for predicting patients’ subtypes. The two circos plots demonstrated how the 30 biomarkers involved in pro-EMT and pro-fibrotic processes. Genes (bottom half of each circle) were linked to the functional annotations (upper half of each circle) that they were involved in.
By doing a literature review on functional records of genes in ImmunCard30, we found that most of the genes were having clear promoting or inhibiting effects on heart failure validated by in-vitro, in-vivo or clinical studies (Supplementary Table 6), and the remaining were found to statistically correlated with onset or progression of heart failure (Supplementary Table 6). Literature review also indicated that functions of 17G-subtype1 and 13G-subtype2 were related to features of their corresponding subtypes. Most genes from 17G-subtype1 directly participated in [for example, TNFRSF17, MEP1B, TNFRSF9, CST6, CANT1, IGF2BP1, IL13, LRP8 (62–71)] or at least were relevant to [for example, PLXNC1, DUOX2, SLITRK3, ST6GAL2, F11 (72–77)] inflammation and adaptive immune responses (Supplementary Table 6), which were major physiological features of the pro-inflammatory subtype. Meanwhile, genes from 13G-subtype2 were related to down-regulation of inflammatory responses [for example, LIMD1, VSIG4, ADAMTS5, VTCN1 (78–81)], blood vessel development [for example, VEGFD, RLN2, KLK12, and ANGPTL3 (82–86)], and maladaptive responses to reduced cardiac functions [for example, VEGFD, ANGPTL3, LIMD1, VSIG4, IDI2 (79, 83–88) (Supplementary Table 6), which were physiological features related to M2 macrophage (13). These records from in-vitro, in-vivo and clinical studies provide further evidence to prove the close relationship between ImmunCard30 and subtypes of immune cell fractions.
In order to evaluate the possibility of applying the ImmunCard30 panel in liquid biopsy, we went through the former literatures related to detection of ImmunCard30 proteins or genes’ expression in blood. According to the result, 27 of the members in ImmunCard30 have clear records indicating correlations between their levels in plasma/blood cells and cardiovascular healthy, while the remaining also have proven correlations with heart failure and being detectable in blood (Supplementary Table 6). Therefore, liquid biopsy based on blood samples would be a potential approach for detecting ImmunCard30. Further studies correlating ventricular immune cell fractions and blood levels of ImmunCard30 would be worth doing to validate the panel.
Overall, we presented a unique perspective based on immune subtyping using immune cell fractions calculated by CIBERSORTx. Subtype affiliation of immune cell fractions was potentially linked to diverse disease progression mechanisms, according to the analyses performed in this study. The accurate stratification of patient subtypes will greatly benefit the development and application of novel personalised therapeutic strategies. The ImmunCard30, which performed well on our discovery and validating sets (n = 280), is a promising biomarker panel that merits further investigation to determine its validity in predicting subtypes for patients and guiding individualised drug therapies.
Discussion
Heart failure is a significant challenge to longevity with complex etiologies, and its progression is regulated by a heterogeneous immune cell fractions and complex feedbacks (1, 18, 26). The maladaptive structural changes in the heart during the progression of heart failure are generalized as cardiac remodeling. Cardiovascular fibrosis is a critical process in cardiac remodeling. Fibrosis was activated as a response to cardiac stress, such as cardiomyocyte death, pressure overload, and myocardial inflammation, to heal dead tissues and improve the physical strength of the myocardium (8). However, in chronic heart failure, fibrotic responses are constantly activated, resulting in excessive extracellular matrix accumulation and pathological fibrosis. The fibrotic tissues would stiffen the ventricle, reducing its compliance and diastolic function. Meanwhile, fibrotic tissues disrupt cardiac electrophysiology and limit nutrient supplies to cardiomyocytes, resulting in abnormal cardiac rhythm and cell death (8).
Because of the central role of fibrosis in the progression of cardiac remodeling, clinical studies have found that fibrosis is strongly associated with the prognosis of heart failure, regardless of the etiology of heart failure (4–6). Therefore, for the treatment of heart failure, it is pivotal to focus on central pathological processes such as cardiac fibrosis (89, 90). However, recent studies of anti-fibrosis therapies targeting fibrotic pathways yielded inconsistent or even worsened outcomes (19, 20), in part attributed to the highly heterogeneous immune cell fractions of fibrotic hearts. Although significant advancements in anti-fibrosis therapies were achieved, the heterogeneity of cell fractions in the cardiac tissue often led to varied treatment responses (13, 91). Precise stratification of patients and individualising therapies considering patients’ heterogeneous would be the future direction (91, 92). Numerous studies have demonstrated that the cardiac immune cell fractions is associated with the mechanism of fibrosis in patients (93, 94). Here, we developed a novel patient stratifying method based on the immune cell fractions and demonstrated the potential of this stratifying method in predicting patients’ fibrotic mechanisms and assisting the design of drug treatments.
We identified two subtypes of immune cell fractions for end-stage heart failure (pro-inflammatory and pro-remodeling subtypes) using ventricular tissue samples from 103 patients with heart failure, with significantly different levels of biomarkers related to adaptive immune responses. In both subtypes, EMT was identified as a potentially central fibrotic process. Endothelial biomarker analyses revealed strong correlations between endothelial cells and EMT. EndMT has been shown in previous studies to be a pro-fibrotic process in disease models (58). Further research on clinical cohorts could help determine whether EndMT was occurring and infer the potential contributions of EndMT to fibrosis in patients with end-stage heart failure. LAFSAA-based evaluation further identified 104 subtype-specific pro-EMT gene sets categorised into 8 functional clusters. We demonstrated that the 11 subtype-specific pro-fibrotic gene sets derived from the 8 functional clusters were related to inflammation, cell-cell adhesion, vasculature development, and responses to environmental stresses. A unique 30-gene biomarker panel was established to accurately stratify patient subtypes in both the discovery set (AUC = 0.954) and the validating set GSE57338 (AUC = 0.803).
Furthermore, subtype-specific pro-fibrotic mechanisms were highly correlated with representative immune cells of both subtypes, indicating casual relationships between the heterogeneous immune cell fractions and fibrotic mechanisms. We demonstrated that multiple processes related to regeneration, blood vessel development, and morphogenesis were correlated with fibrosis in the pro-remodeling subtype [for example, varicosity, smooth muscle cell matrix adhesion (Figures 4A, B)], while M2 macrophages (a distinct phenotype of the pro-remodeling subtype) was highly related to these processes (13, 61). In contrast, inflammation and adaptive immune responses, which were primarily modulated by T cells, B cells, and monocytes, were specifically correlated with fibrosis in the pro-inflammatory subtype. For example, defense response to virus, a typical adaptive immune response, was a subtype-specific fibrotic process in the pro-inflammatory subtype (Figure 4B). IL-13 and IL-4 secreted by CD4+ T cells (a distinct phenotype of the pro-inflammatory subtype) activate Janus kinase activity (95), and activation of Janus kinase was identified as a subtype-specific fibrotic process in the pro-inflammatory subtype (Figure 4B). These interplays between fibrosis and immune cells were also observed in early-stage heart failure (96–98). Therefore, although our current results were obtained based on end-stage patients, the panel should also be applicable to early-stage patients.
We envision that the 30-gene biomarker panel derived from ventricular tissue samples could provide a key reference for patient stratification, facilitating personalised therapeutic strategies based on drugs targeting fibrotic pathways (especially those related to subtype-specific pro-fibrotic gene sets, such as endothelin inhibitors, anti-TGFβ drugs, and anti-inflammation drugs), and reduce the adverse symptoms caused by drug treatments (such as inflammatory response (21)). Although ventricular tissue samples can be obtained with endomyocardial biopsy (99), the procedure increases mortality risk, cost and inconveniences the patient. Routine screening and assessment of patient subtypes would be more feasible through minimally invasive methods such as cardiovascular magnetic resonance (CMR), which is non-invasive and is the primary technology for stratifying patients with cardiac fibrosis (19, 100). However, CMR cannot provide information related to the intrinsic immunity of fibrotic lesions and is targeted to monitor fibrosis rather than guide individualised drug treatments targeting the fibrotic pathways. Here, our novel strategy for assessing subtypes of immune cell fractions is highly relevant to disease progression mechanisms. The ImmunCard30 is the first diagnostic biomarker panel reported to potentially reflect the heart’s immune cell fractions, assisting in assessing patients’ fibrotic mechanisms in the heart. The ImmunCard30 biomarker panel’s correlation with patients’ pathological conditions is expected to have significant implications for personalised medicine and, as such, warrants further clinical studies to validate, particularly validation using minimally invasive assessment methods such as liquid biopsy.
Clinical sampling of ventricular tissues is challenging due to the high risk of sampling operation. Immuncard30 is based on data from failing heart samples (n = 208) from various sources. However, there could still be bias in evaluating subtype-association of fibrotic mechanisms and the ImmunCard30 biomarker panel’s classifying efficiency. Our study provided new insights into the pathological correlations of cardiac immune cell fractions, which will be of great interest in developing personalised medicine. Further studies focusing on the in vitro and in vivo validation of the biomarker panel to validate drug responses and diagnostic efficiency are in progress to promote actual clinical utility. The progress made by this study support effort for the development of anti-fibrosis precision medicine for heart failure.
Data availability statement
The original contributions presented in the study are included in the article/Supplementary Material. Further inquiries can be directed to the corresponding author.
Author contributions
SZ conceived and designed the study, performed the study, analysed the data, and wrote the paper. BK conceived and designed the study, and wrote the paper. All authors contributed to the article and approved the submitted version.
Funding
This work was supported in part by InnoHK Project on [Project 1-2 - Design of DDS system structure and selection of potential drugs for acute CVDs] and [Project 1-5 Multi-modal spectroscopy (MMS) & biosensor platforms for monitoring CVDs] at Hong Kong Centre for Cerebro-cardiovascular Health Engineering (COCHE), in part by City University of Hong Kong (9610430, 7020002, 7005464, 7005208, 9667220), which is funded by the Research Grants Council (RGC), in part by Pneumoconiosis Compensation Fund Board (9211276), and in part by Research Grants Council (RGC 9048206). The views expressed are those of the authors and not necessarily those of InnoHK - ITC, City University of Hong Kong, Pneumoconiosis Compensation Fund Board, or Research Grants Council.
Conflict of interest
The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.
Publisher’s note
All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.
Supplementary material
The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fimmu.2023.1053793/full#supplementary-material
References
1. Follath F, Cleland JGF, Klein W, Murphy R. Etiology and response to drug treatment in heart failure. J Am Coll Cardiol (1998) 32(5):1167–72. doi: 10.1016/S0735-1097(98)00400-8
2. Schultheiss H-P, Fairweather D, Caforio ALP, Escher F, Hershberger RE, Lipshultz SE, et al. Dilated cardiomyopathy. Nat Rev Dis Primers (2019) 5(1):1–19. doi: 10.1038/s41572-019-0084-1
3. Bacmeister L, Schwarzl M, Warnke S, Stoffers B, Blankenberg S, Westermann D, et al. Inflammation and fibrosis in murine models of heart failure. Basic Res Cardiol (2019) 114(3):19–35. doi: 10.1007/s00395-019-0722-5
4. Ahmad T, Fiuzat M, Neely B, Neely ML, Pencina MJ, Kraus WE, et al. Biomarkers of myocardial stress and fibrosis as predictors of mode of death in patients with chronic heart failure. JACC: Heart Failure (2014) 2(3):260–8. doi: 10.1016/j.jchf.2013.12.004
5. Mandawat A, Chattranukulchai P, Mandawat A, Blood AJ, Ambati S, Hayes B, et al. Progression of myocardial fibrosis in nonischemic dcm and association with mortality and heart failure outcomes. JACC Cardiovasc Imaging (2021) 14(7):1338–50. doi: 10.1016/j.jcmg.2020.11.006
6. Gulati A, Japp AG, Raza S, Halliday BP, Jones DA, Newsome S, et al. Absence of myocardial fibrosis predicts favorable long-term survival in new-onset heart failure. Circulation: Cardiovasc Imaging (2018) 11(9):e007722. doi: 10.1161/CIRCIMAGING.118.007722
7. Ichiki T, Schirger JA, Huntley BK, Brozovich FV, Maleszewski JJ, Sandberg SM, et al. Cardiac fibrosis in end-stage human heart failure and the cardiac natriuretic peptide guanylyl cyclase system: Regulation and therapeutic implications. J Mol Cell Cardiol (2014) 75:199–205. doi: 10.1016/j.yjmcc.2014.08.001
8. Piek A, de Boer RA, Silljé HHW. The fibrosis-cell death axis in heart failure. Heart Failure Rev (2016) 21(2):199–211. doi: 10.1007/s10741-016-9536-9
9. Rossignol P, Hernandez AF, Solomon SD, Zannad F. Heart failure drug treatment. Lancet (2019) 393(10175):1034–44. doi: 10.1016/S0140-6736(18)31808-7
10. Yancy CW, Jessup M, Bozkurt B, Butler J, Casey JDE, Colvin MM, et al. 2017 Acc/Aha/Hfsa focused update of the 2013 Accf/Aha guideline for the management of heart failure: A report of the American college of Cardiology/American heart association task force on clinical practice guidelines and the heart failure society of America. J Am Coll Cardiol (2017) 70(6):776–803. doi: 10.1016/j.jacc.2017.04.025
11. Liu M, Lopez de Juan Abad B, Cheng K. Cardiac fibrosis: Myofibroblast-mediated pathological regulation and drug delivery strategies. Adv Drug Delivery Rev (2021) 173:504–19. doi: 10.1016/j.addr.2021.03.021
12. Rurik JG, Aghajanian H, Epstein JA. Immune cells and immunotherapy for cardiac injury and repair. Circ Res (2021) 128(11):1766–79. doi: 10.1161/Circresaha.121.318005
13. Sun K, Li YY, Jin J. A double-edged sword of immuno-microenvironment in cardiac homeostasis and injury repair. Signal Transduct Tar (2021) 6(1):79–. doi: 10.1038/s41392-020-00455-6
14. Blanton RM, Carrillo-Salinas FJ, Alcaide P. T-Cell recruitment to the heart: Friendly guests or unwelcome visitors? Am J Physiol Heart Circ Physiol (2019) 317(1):H124–H40. doi: 10.1152/ajpheart.00028.2019
15. Hou XZ, Chen GB, Bracamonte-Baran W, Choi HS, Diny NL, Sung J, et al. The cardiac microenvironment instructs divergent monocyte fates and functions in myocarditis. Cell Rep (2019) 28(1):172–89.e7. doi: 10.1016/j.celrep.2019.06.007
16. Lis-López L, Bauset C, Seco-Cervera M, Cosín-Roger J. Is the macrophage phenotype determinant for fibrosis development? Biomedicines (2021) 9(12):1747. doi: 10.3390/biomedicines9121747
17. Okyere AD, Tilley DG. Leukocyte-dependent regulation of cardiac fibrosis. Front Physiol (2020) 11:301. doi: 10.3389/fphys.2020.00301
18. Aujla PK, Kassiri Z. Diverse origins and activation of fibroblasts in cardiac fibrosis. Cell signalling (2021) 78:109869. doi: 10.1016/j.cellsig.2020.109869
19. Webber M, Jackson SP, Moon JC, Captur G. Myocardial fibrosis in heart failure: Anti-fibrotic therapies and the role of cardiovascular magnetic resonance in drug trials. Cardiol Ther (2020) 9(2):363–76. doi: 10.1007/s40119-020-00199-y
20. Murphy SP, Kakkar R, McCarthy CP, Januzzi JL. Inflammation in heart failure jacc state-of-the-Art review. J Am Coll Cardiol (2020) 75(11):1324–40. doi: 10.1016/j.jacc.2020.01.014
21. Engebretsen KVT, Skardal K, Bjornstad S, Marstein HS, Skrbic B, Sjaastad I, et al. Attenuated development of cardiac fibrosis in left ventricular pressure overload by Sm16, an orally active inhibitor of Alk5. J Mol Cell Cardiol (2014) 76:148–57. doi: 10.1016/j.yjmcc.2014.08.008
22. Frangogiannis NG. Cardiac fibrosis: Cell biological mechanisms, molecular pathways and therapeutic opportunities. Mol aspects Med (2019) 65:70–99. doi: 10.1016/j.mam.2018.07.001
23. Wei L. Immunological aspect of cardiac remodeling: T lymphocyte subsets in inflammation-mediated cardiac fibrosis. Exp Mol Pathol (2011) 90(1):74–8. doi: 10.1016/j.yexmp.2010.10.004
24. García-Rivas G, Castillo EC, Gonzalez-Gil AM, Maravillas-Montero JL, Brunck M, Torres-Quintanilla A, et al. The role of b cells in heart failure and implications for future immunomodulatory treatment strategies. ESC Heart Fail (2020) 7(4):1387–99. doi: 10.1002/ehf2.12744
25. Weldy CS, Ashley EA. Towards precision medicine in heart failure. Nat Rev Cardiol (2021) 18(11):745–62. doi: 10.1038/s41569-021-00566-9
26. Leopold JA, Maron BA, Loscalzo J. The application of big data to cardiovascular disease: Paths to precision medicine. J Clin Invest (2020) 130(1):29–38. doi: 10.1172/Jci129203
27. Yamaguchi T, Sumida TS, Nomura S, Satoh M, Higo T, Ito M, et al. Cardiac dopamine D1 receptor triggers ventricular arrhythmia in chronic heart failure. Nat Commun (2020) 11(1):4364–. doi: 10.1038/s41467-020-18128-x
28. Hua X, Wang YY, Jia P, Xiong Q, Hu Y, Chang Y, et al. Multi-level transcriptome sequencing identifies Col1a1 as a candidate marker in human heart failure progression. BMC Med (2020) 18(1):2–. doi: 10.1186/s12916-019-1469-4
29. Yang KC, Yamada KA, Patel AY, Topkara VK, George I, Cheema FH, et al. Deep rna sequencing reveals dynamic regulation of myocardial noncoding rnas in failing human heart and remodeling with mechanical circulatory support. Circulation (2014) 129(9):1009–21. doi: 10.1161/CIRCULATIONAHA.113.003863
30. Liu YC, Morley M, Brandimarto J, Hannenhalli S, Hu Y, Ashley EA, et al. Rna-seq identifies novel myocardial gene expression signatures of heart failure. Genomics (2015) 105(2):83–9. doi: 10.1016/j.ygeno.2014.12.002
31. McDonough JE, Ahangari F, Li Q, Jain S, Verleden SE, Herazo-Maya J, et al. Transcriptional regulatory model of fibrosis progression in the human lung. JCI Insight (2019) 4(22):2379–3708. doi: 10.1172/jci.insight.131597
32. Robinson MD, McCarthy DJ, Smyth GK. Edger: A bioconductor package for differential expression analysis of digital gene expression data. Bioinformatics (2010) 26(1):139–40. doi: 10.1093/bioinformatics/btp616
33. Liberzon A, Birger C, Thorvaldsdóttir H, Ghandi M, Mesirov JP, Tamayo P. The molecular signatures database (Msigdb) hallmark gene set collection. Cell Syst (2015) 1(6):417–25. doi: 10.1016/j.cels.2015.12.004
34. Subramanian A, Tamayo P, Mootha VK, Mukherjee S, Ebert BL, Gillette MA, et al. Gene set enrichment analysis: A knowledge-based approach for interpreting genome-wide expression profiles. P Natl Acad Sci USA (2005) 102(43):15545–50. doi: 10.1073/pnas.0506580102
35. Uhlen M, Karlsson MJ, Hober A, Svensson AS, Scheffel J, Kotol D, et al. The human secretome. Sci Signal (2019) 12(609):eaaz0274. doi: 10.1126/scisignal.aaz0274
36. Newman AM, Liu CL, Green MR, Gentles AJ, Feng W, Xu Y, et al. Robust enumeration of cell subsets from tissue expression profiles. Nat Methods (2015) 12(5):453–7. doi: 10.1038/nmeth.3337
37. Newman AM, Steen CB, Liu CL, Gentles AJ, Chaudhuri AA, Scherer F, et al. Determining cell type abundance and expression from bulk tissues with digital cytometry. Nat Biotechnol (2019) 37(7):773–82. doi: 10.1038/s41587-019-0114-2
38. Rousseeuw PJ. Silhouettes - a graphical aid to the interpretation and validation of cluster-analysis. J Comput Appl Math (1987) 20:53–65. doi: 10.1016/0377-0427(87)90125-7
39. Kassambara A, Mundt F. Extract and visualize the results of multivariate data analyses [R package factoextra version 1.0.7]. (2020). Available at: https://CRAN.R-project.org/package=factoextra
40. Van der Maaten L, Hinton G. Visualizing data using T-sne. J Mach Learn Res (2008) 9(11):2579–605.
41. Kuhn M. Building predictive models in r using the caret package. J Stat software (2008) 28(5):1–26. doi: 10.18637/jss.v028.i05
42. Love MI, Huber W, Anders S. Moderated estimation of fold change and dispersion for rna-seq data with Deseq2. Genome Biol (2014) 15(12):550–. doi: 10.1186/s13059-014-0550-8
43. Zhou YY, Zhou B, Pache L, Chang M, Khodabakhshi AH, Tanaseichuk O, et al. Metascape provides a biologist-oriented resource for the analysis of systems-level datasets. Nat Commun (2019) 10:1523. doi: 10.1038/s41467-019-09234-6
44. Foroutan M, Bhuva DD, Lyu R, Horan K, Cursons J, Davis MJ. Single sample scoring of molecular phenotypes. BMC Bioinf (2018) 19:404. doi: 10.1186/s12859-018-2435-4
45. Friedman J, Hastie T, Tibshirani R. Regularization paths for generalized linear models Via coordinate descent. J Stat Software (2010) 33(1):1–22. doi: 10.18637/jss.v033.i01
46. Robin X, Turck N, Hainard A, Tiberti N, Lisacek F, Sanchez JC, et al. Proc: An open-source package for r and s plus to analyze and compare roc curves. BMC Bioinf (2011) 12:77–. doi: 10.1186/1471-2105-12-77
47. Chi X, Sartor MA, Lee S, Anurag M, Patil S, Hall P, et al. Universal concept signature analysis: Genome-wide quantification of new biological and pathological functions of genes and pathways. Brief Bioinform (2020) 21:1717–32. doi: 10.1093/bib/bbz093
48. Li L, Zhao Q, Kong W. Extracellular matrix remodeling and cardiac fibrosis. Matrix Biol (2018) 68-69:490–506. doi: 10.1016/j.matbio.2018.01.013
49. Heymans S, Schroen B, Vermeersch P, Milting H, Gao F, Kassner A, et al. Increased cardiac expression of tissue inhibitor of metalloproteinase-1 and tissue inhibitor of metalloproteinase-2 is related to cardiac fibrosis and dysfunction in the chronic pressure-overloaded human heart. Circulation (2005) 112(8):1136–44. doi: 10.1161/circulationaha.104.516963
50. Chothani S, Schafer S, Adami E, Viswanathan S, Widjaja AA, Langley SR, et al. Widespread translational control of fibrosis in the human heart by rna-binding proteins. Circulation (2019) 140(11):937–51. doi: 10.1161/CIRCULATIONAHA.119.039596
51. Kaur H, Takefuji M, Ngai CY, Carvalho J, Bayer J, Wietelmann A, et al. Targeted ablation of periostin-expressing activated fibroblasts prevents adverse cardiac remodeling in mice. Circ Res (2016) 118(12):1906–17. doi: 10.1161/circresaha.116.308643
52. Glezeva N, Horgan S, Baugh JA. Monocyte and macrophage subsets along the continuum to heart failure: Misguided heroes or targetable villains? J Mol Cell Cardiol (2015) 89:136–45. doi: 10.1016/j.yjmcc.2015.10.029
53. Roberts CA, Dickinson AK, Taams LS. The interplay between Monocytes/Macrophages and Cd4(+) T cell subsets in rheumatoid arthritis. Front Immunol (2015) 6:571. doi: 10.3389/fimmu.2015.00571
54. Kovacic JC, Mercader N, Torres M, Boehm M, Fuster V. Epithelial-to-Mesenchymal and endothelial-to-Mesenchymal transition: From cardiovascular development to disease. Circulation (2012) 125(14):1795–808. doi: 10.1161/circulationaha.111.040352
55. Shu DY, Butcher E, Saint-Geniez M. Emt and endmt: Emerging roles in age-related macular degeneration. Int J Mol Sci (2020) 21(12):4271–. doi: 10.3390/ijms21124271
56. Yang J, Antin P, Berx G, Blanpain C, Brabletz T, Bronner M, et al. Guidelines and definitions for research on epithelial–mesenchymal transition. Nat Rev Mol Cell Biol (2020) 21(6):341–52. doi: 10.1038/s41580-020-0237-9
57. Piera-Velazquez S, Jimenez SA. Endothelial to mesenchymal transition: Role in physiology and in the pathogenesis of human diseases. Physiol Rev (2019) 99(2):1281–324. doi: 10.1152/physrev.00021.2018
58. Zeisberg EM, Tarnavski O, Zeisberg M, Dorfman AL, McMullen JR, Gustafsson E, et al. Endothelial-to-Mesenchymal transition contributes to cardiac fibrosis. Nat Med (2007) 13(8):952–61. doi: 10.1038/nm1613
59. Meyer-Schaller N, Cardner M, Diepenbruck M, Saxena M, Tiede S, Luond F, et al. A hierarchical regulatory landscape during the multiple stages of emt. Dev Cell (2019) 48(4):539–53 e6. doi: 10.1016/j.devcel.2018.12.023
60. Villar-Fincheira P, Sanhueza-Olivares F, Norambuena-Soto I, Cancino-Arenas N, Hernandez-Vargas F, Troncoso R, et al. Role of interleukin-6 in vascular health and disease. Front Mol Biosci (2021) 8:641734. doi: 10.3389/fmolb.2021.641734
61. Hong H, Tian XY. The role of macrophages in vascular repair and regeneration after ischemic injury. Int J Mol Sci (2020) 21(17):6328–. doi: 10.3390/ijms21176328
62. Saunders CN, Chattopadhyay S, Huhn S, Weinhold N, Hoffmann P, Nöthen MM, et al. Search for Al amyloidosis risk factors using mendelian randomization. Blood Adv (2021) 5(13):2725–31. doi: 10.1182/bloodadvances.2021004423
63. Frigyesi I, Adolfsson J, Ali M, Kronborg Christophersen M, Johnsson E, Turesson I, et al. Robust isolation of malignant plasma cells in multiple myeloma. Blood (2014) 123(9):1336–40. doi: 10.1182/blood-2013-09-529800
64. Sun Q, Jin HJ, Bond JS. Disruption of the meprin A and B genes in mice alters homeostasis of monocytes and natural killer cells. Exp Hematol (2009) 37(3):346–56. doi: 10.1016/j.exphem.2008.10.016
65. Li Y-J, Fan Y-H, Tang J, Li J-B, Yu C-H. Meprin-B regulates production of pro-inflammatory factors Via a disintegrin and metalloproteinase-10 (Adam-10) dependent pathway in macrophages. Int Immunopharmacol (2014) 18(1):77–84. doi: 10.1016/j.intimp.2013.11.004
66. Olofsson PS. Targeting T cell costimulation to prevent atherothrombosis. Circ Res (2012) 110(6):800–1. doi: 10.1161/CIRCRESAHA.112.265108
67. Zannad F, Ferreira JP, Butler J, Filippatos G, Januzzi JL, Sumin M, et al. Effect of empagliflozin on circulating proteomics in heart failure: Mechanistic insights into the emperor programme. Eur Heart J (2022) 43(48):4991–5002. doi: 10.1093/eurheartj/ehac495
68. Gao F, Hu X, Liu W, Wu H, Mu Y, Zhao Y. Calcium-activated nucleotides 1 (Cant1)-driven nuclear factor-K-Gene binding (Nf-ĸb) signaling pathway facilitates the lung cancer progression. Bioengineered (2022) 13(2):3183–93. doi: 10.1080/21655979.2021.2003131
69. Xie J, Li Q, Zhu X-H, Gao Y, Zhao W-H. Igf2bp1 promotes lps-induced nfκb activation and pro-inflammatory cytokines production in human macrophages and monocytes. Biochem Biophys Res Commun (2019) 513(4):820–6. doi: 10.1016/j.bbrc.2019.03.206
70. Qian N, Gao Y, Wang J, Wang Y. Emerging role of interleukin-13 in cardiovascular diseases: A ray of hope. J Cell Mol Med (2021) 25(12):5351–7. doi: 10.1111/jcmm.16566
71. Wei J, Cui J. Mir-1322 protects against the myocardial ischemia Via Lrp8/Pi3k/Akt pathway. Biochem Biophys Res Commun (2023) 638:120–6. doi: 10.1016/j.bbrc.2022.10.101
72. König K, Marth L, Roissant J, Granja T, Jennewein C, Devanathan V, et al. The plexin C1 receptor promotes acute inflammation. Eur J Immunol (2014) 44(9):2648–58. doi: 10.1002/eji.201343968
73. Zhang H, Bredewold EOW, Vreeken D, Duijs JMGJ, de Boer HC, Kraaijeveld AO, et al. Prediction power on cardiovascular disease of neuroimmune guidance cues expression by peripheral blood monocytes determined by machine-learning methods. Int J Mol Sci (2020) 21(17):6364. doi: 10.3390/ijms21176364
74. Milanesi E, Manda G, Dobre M, Codrici E, Neagoe IV, Popescu BO, et al. Distinctive under-expression profile of inflammatory and redox genes in the blood of elderly patients with cardiovascular disease. J Inflammation Res (2021) 14:429–42. doi: 10.2147/jir.S280328
75. Afshari P, Yao WD, Middleton FA. Reduced Slc1a1 expression is associated with neuroinflammation and impaired sensorimotor gating and cognitive performance in mice: Implications for schizophrenia. PloS One (2017) 12(9):e0183854. doi: 10.1371/journal.pone.0183854
76. Verdonschot JAJ, Wang P, Van Bilsen M, Hazebroek MR, Merken JJ, Vanhoutte EK, et al. Metabolic profiling associates with disease severity in nonischemic dilated cardiomyopathy. J Cardiac Failure (2020) 26(3):212–22. doi: 10.1016/j.cardfail.2019.09.004
77. Cao Y, Wang Y, Zhou Z, Pan C, Jiang L, Zhou Z, et al. Liver-heart cross-talk mediated by coagulation factor xi protects against heart failure. Science (2022) 377(6613):1399–406. doi: 10.1126/science.abn0910
78. Garcia-Campos A, Correia CN, Naranjo-Lucena A, Garza-Cuartero L, Farries G, Browne JA, et al. Fasciola hepatica infection in cattle: Analyzing responses of peripheral blood mononuclear cells (Pbmc) using a transcriptomics approach. Front Immunol (2019) 10:2081. doi: 10.3389/fimmu.2019.02081
79. di Salvo TG, Yang K-C, Brittain E, Absi T, Maltais S, Hemnes A. Right ventricular myocardial biomarkers in human heart failure. J Cardiac Failure (2015) 21(5):398–411. doi: 10.1016/j.cardfail.2015.02.005
80. Barallobre-Barreiro J, Radovits T, Fava M, Mayr U, Lin W-Y, Ermolaeva E, et al. Extracellular matrix in heart failure: Role of Adamts5 in proteoglycan remodeling. Circulation (2021) 144(25):2021–34. doi: 10.1161/CIRCULATIONAHA.121.055732
81. Vaishnav J, Khan F, Yadav M, Parmar N, Buch H, Jadeja SD, et al. V-Set domain containing T-cell activation inhibitor-1 (Vtcn1): A potential target for the treatment of autoimmune diseases. Immunobiology (2022) 227(6):152274. doi: 10.1016/j.imbio.2022.152274
82. Kryza T, Parent C, Pardessus J, Petit A, Burlaud-Gaillard J, Reverdiau P, et al. Human kallikrein-related peptidase 12 stimulates endothelial cell migration by remodeling the fibronectin matrix. Sci Rep (2018) 8(1):6331. doi: 10.1038/s41598-018-24576-9
83. Zhou Y, Zhu X, Cui H, Shi J, Yuan G, Shi S, et al. The role of the vegf family in coronary heart disease. Front Cardiovasc Med (2021) 8:738325. doi: 10.3389/fcvm.2021.738325
84. Iguchi M, Wada H, Shinozaki T, Suzuki M, Ajiro Y, Matsuda M, et al. Distinct association of vegf-c and vegf-d with prognosis in patients with chronic heart failure: The prehosp-chf study. Eur Heart J (2021) 42(Supplement_1). doi: 10.1093/eurheartj/ehab724.0868
85. Camenisch G, Pisabarro MT, Sherman D, Kowalski J, Nagel M, Hass P, et al. Angptl3 stimulates endothelial cell adhesion and migration Via integrin migration Via integrin Avβ3 and induces blood vessel formation in vivo. J Biol Chem (2002) 277(19):17281–90. doi: 10.1074/jbc.M109768200
86. Ponikowski P, Metra M, Teerlink JR, Unemori E, Felker GM, Voors AA, et al. Design of the relaxin in acute heart failure study. Am Heart J (2012) 163(2):149–55.e1. doi: 10.1016/j.ahj.2011.10.009
87. Landry NM, Rattan SG, Filomeno KL, Meier TW, Meier SC, Foran SJ, et al. Ski activates the hippo pathway Via Limd1 to inhibit cardiac fibroblast activation. Basic Res Cardiol (2021) 116(1):25. doi: 10.1007/s00395-021-00865-9
88. Luo Y, Mao C, Yang Y, Wang F, Ahmad FS, Arnett D, et al. Integrating hypertension phenotype and genotype with hybrid non-negative matrix factorization. Bioinformatics (2018) 35(8):1395–403. doi: 10.1093/bioinformatics/bty804
89. de Boer RA, De Keulenaer G, Bauersachs J, Brutsaert D, Cleland JG, Diez J, et al. Towards better definition, quantification and treatment of fibrosis in heart failure. a scientific roadmap by the committee of translational research of the heart failure association (Hfa) of the European society of cardiology. Eur J Heart Fail (2019) 21(3):272–85. doi: 10.1002/ejhf.1406
90. Travers JG, Tharp CA, Rubino M, McKinsey TA. Therapeutic targets for cardiac fibrosis: From old school to next-gen. J Clin Invest (2022) 132(5):1–11. doi: 10.1172/JCI148554
91. Henderson NC, Rieder F, Wynn TA. Fibrosis: From mechanisms to medicines. Nature (2020) 587(7835):555–66. doi: 10.1038/s41586-020-2938-9
92. Prabhu SD, Frangogiannis NG. The biological basis for cardiac repair after myocardial infarction from inflammation to fibrosis. Circ Res (2016) 119(1):91–112. doi: 10.1161/Circresaha.116.303577
93. Liu X, Shi GP, Guo JL. Innate immune cells in pressure overload-induced cardiac hypertrophy and remodeling. Front Cell Dev Biol (2021) 9:659666. doi: 10.3389/fcell.2021.659666
94. Theall B, Alcaide P. The heart under pressure: Immune cells in fibrotic remodeling. Curr Opin Physiol (2022) 25:100484. doi: 10.1016/j.cophys.2022.100484
95. Nguyen JK, Austin E, Huang A, Mamalis A, Jagdeo J. The il-4/Il-13 axis in skin fibrosis and scarring: Mechanistic concepts and therapeutic targets. Arch Dermatol Res (2020) 312(2):81–92. doi: 10.1007/s00403-019-01972-3
96. Shiraishi M, Shintani Y, Shintani Y, Ishida H, Saba R, Yamaguchi A, et al. Alternatively activated macrophages determine repair of the infarcted adult murine heart. J Clin Invest (2016) 126(6):2151–66. doi: 10.1172/Jci85782
97. Stawski L, Haines P, Fine A, Rudnicka L, Trojanowska M. Mmp-12 deficiency attenuates angiotensin ii-induced vascular injury, M2 macrophage accumulation, and skin and heart fibrosis. PloS One (2014) 9(10):e109763. doi: 10.1371/journal.pone.0109763
98. Laroumanie F, Douin-Echinard V, Pozzo J, Lairez O, Tortosa F, Vinel C, et al. Cd4(+) T cells promote the transition from hypertrophy to heart failure during chronic pressure overload. Circulation (2014) 129(21):2111–24. doi: 10.1161/Circulationaha.113.007101
99. Xiu MX, Liu YM, Wang WJ. Investigation of hub genes and immune status in heart transplant rejection using endomyocardial biopsies. J Cell Mol Med (2021) 25(2):763–73. doi: 10.1111/jcmm.16127
100. Grani C, Eichhorn C, Biere L, Kaneko K, Murthy VL, Agarwal V, et al. Comparison of myocardial fibrosis quantification methods by cardiovascular magnetic resonance imaging for risk stratification of patients with suspected myocarditis. J Cardiovasc Magn R (2019) 21:14–. doi: 10.1186/s12968-019-0520-0
Keywords: heart failure, cardiac fibrosis, immune cell fractions, subtyping, functional gene sets
Citation: Zou S and Khoo BL (2023) Subtyping based on immune cell fractions reveal heterogeneity of cardiac fibrosis in end-stage heart failure. Front. Immunol. 14:1053793. doi: 10.3389/fimmu.2023.1053793
Received: 26 September 2022; Accepted: 02 February 2023;
Published: 15 February 2023.
Edited by:
Chris Wincup, King’s College Hospital NHS Foundation Trust, United KingdomReviewed by:
Qiyuan Shan, Zhejiang Chinese Medical University, ChinaDaniel Abebayehu, University of Virginia, United States
Copyright © 2023 Zou and Khoo. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.
*Correspondence: Bee Luan Khoo, YmxraG9vQGNpdHl1LmVkdS5oaw==