Characterization of the transcriptional response of Candida parapsilosis to the antifungal peptide MAF-1A

Candida parapsilosis is a major fungal pathogen that leads to sepsis. New and more effective antifungal agents are required due to the emergence of resistant fungal strains. MAF-1A is a cationic antifungal peptide isolated from Musca domestica that is effective against a variety of Candida species. However, the mechanism(s) of its antifungal activity remains undefined. Here, we used RNA-seq to identify differentially expressed genes (DEGs) in Candida parapsilosis following MAF-1A exposure. The early (6 h) response included 1,122 upregulated and 1,065 downregulated genes. Late (18 h) responses were associated with the increased expression of 101 genes and the decreased expression of 151 genes. Upon MAF-1A treatment for 18 h, 42 genes were upregulated and 25 genes were downregulated. KEGG enrichment showed that the DEGs in response to MAF-1A were mainly involved in amino acid synthesis and metabolism, oxidative phosphorylation, sterol synthesis, and apoptosis. These results indicate that MAF-1A exerts antifungal activity through interference with Candida parapsilosis cell membrane integrity and organelle function. This provides new insight into the interaction between Candida parapsilosis and this antimicrobial peptide and serves as a reference for future Candida parapsilosis therapies.


INTRODUCTION
Immunosuppressed patients are at a high risk of hospital-acquired fungal infections. Candida albicans (C. albicans) is the most common pathogen of Candida species, its dominance has decreased as the incidence of non-albicans Candida (NAC) species have increased (Vieira de Melo et al., 2019). Over the last two decades, epidemiological studies of fungal pathogens have shown that NAC has surpassed C. albicans as the most prevalent cause of invasive Candida (Sular et al., 2018). New anti-NAC treatment regimens are therefore urgently required.
In recent years, the development of high throughput sequencing technologies has facilitated research on both antimicrobial drug function and drug-resistance. For example, HAC1 (CPAR2_103720) is a key mediator of endoplasmic reticulum stress in C. parapsilosis identified through RNA-seq (Iracane et al., 2018). In our previous studies, we showed that MAF-1A inhibits C. albicans through its effects on the cell wall, plasma membrane, protein synthesis, and energy metabolism (Wang et al., 2017). However, the mechanism(s) through which C. albicans responds to MAF-1A were not fully defined. Here, we have expanded our knowledge on how MAF-1A acts on C. parapsilosis and investigated differences in the responses of C. albicans and C. parapsilosis to MAF-1A treatment. RNA-seq was used to investigate changes in gene expression at early (6 h) and late (18 h) time points, according to time-kill curves of C. parapsilosis growth.

Strains and growth conditions
Transcriptional profiling was performed on the C. parapsilosis reference strain ATCC22019. The strain was preserved in goat blood and stored at −80 • C. C. parapsilosis was streaked on Sabouraud Dextrose Agar (SDA) plates (Sangon, Shanghai, China) at 35 • C as described by Lis et al. (2010). MAF-1A treatments were performed in Sabouraud Dextrose Broth (SDB) (Sangon, Shanghai, China).

Peptide synthesis
MAF-1A was synthesized by Sangon Biotech (Shanghai, Shanghai, China) as a linear peptide of 26 amino acids: KKFKETADKLIESAKQQLESLAKEMK. Analytical high-performance liquid chromatography (HPLC) was used to confirm purity ≥ 95%. The peptide was dissolved in sterile ultrapure water at 5 mg/mL and stored at −20 • C.

Minimum inhibitory concentration (MIC) and time-kill curves
Antifungal assays were performed as per the requirements of the Clinical and Laboratory Standards Institute (CLSI) M27-A3. Briefly, cultures were grown for 24 h at 35 • C and resuspended in SDB. Concentrations were adjusted to approximately 0.5 × 10 3 -2.5 × 10 3 CFU/mL and 100 µl of the suspension was added to each well of 96-well polypropylene microplates (NEST, Wuxi, China). MAF-1A was added at concentration ranging from 0.1 mg/mL to 1.2 mg/mL. All experiments were performed in triplicates. After incubation at 35 • C for 24 h, absorbances were measured at 492 nm on a Microplate Reader (BioTek Synergy H1, Vermont, USA). MIC was defined as the lowest drug concentration showing 80% growth inhibition compared to the drug-free controls. The following formulas were used (Li et al., 2008): (2) Percentage Inhibition of Fungal Growth = 1 − Percentage Fungal Growth. Time-kill curves were performed according to the literature (Li et al., 2008;Sun et al., 2008). C. parapsilosis suspensions were mixed with MAF-1A (MIC) in triplicate and cultured at 35 • C. Aliquots of 100 µl were removed from each test solution at predetermined time points (0, 2, 4, 6, 8, 10, 12, 14, 16, 18, 20, 22, and 24 h). Dilutions were produced (1:100) and streaked in triplicate onto SDA agar plates for colony counts after incubation at 35 • C for 24 h. Sterile ultrapure water was used as a control.

Transcriptome sequencing
C. parapsilosis was inoculated into SDB medium (Sangon, Shanghai, China) at 35 • C for 24 h. C. parapsilosis was treated with MAF-1A at MIC for 6 h (CPAS) and 18 h (CPBS), before RNA extraction. Untreated cultures served as controls (6 h, CPAC; 18 h, CPBC). Total RNA was extracted using RNAiso Plus (Takara, Dalian, China) according to the manufacturer's instructions. RNA concentration and quality were determined on a NanoDrop 2000 (Thermo Fisher Scientific, Wilmington, DE, USA) and Agilent 2100 bioanalyzer (Agilent Technologies, CA, USA). Libraries were prepared using NEBNext R UltraTM RNA Library Prep Kit (NEB, USA) as per the manufacturer's recommendations. Purified libraries were quantified on an Agilent 2100 bioanalyzer Effective concentrations were determined through qRT-PCR analysis. Libraries were prepared and sequenced using a Novoseq sequencer (Illumina, USA) to produce 150 bp paired-end reads.

Differential expression analysis
Raw reads were filtered to obtain high-quality clean reads for subsequent analysis. All reads were mapped to the reference genome of C. parapsilosis (assembly ASM18276v2) from the National Center for Biotechnology Information (NCBI) using HISAT2 v2.0.5 (Kim, Langmead & Salzberg, 2015). Differential expression analysis between the conditions was assessed using the Bioconductor software package DESeq2 in R 1.16.1 (Love, Huber & Anders, 2014). Relative gene expression was assessed using FPKM (Fragments Per Kilobase of transcript sequence per Millions of base pairs sequenced) and compared using log2 FC. P-values were adjusted to generate false discovery rates (padj) as described by Benjamini-Yekutieli, assigning the significance threshold for DEGs as padj < 0.05 (Benjamini et al., 2001;Mortazavi et al., 2008).

Enrichment and interaction network analysis of the differentially expressed genes
To further understand the functions of the DEGs, Gene Ontology (GO) enrichment was performed using the Bioconductor software clusterProfiler 3.4.4 in R package (Yu et al., 2012). Statistical enrichment of the DEGs was also performed in the Kyoto Encyclopedia of Genes and Genomes (KEGG) for pathway enrichment (Kanehisa et al., 2019;Ogata et al., 1999). PPI analysis of the DEGs was performed based on the STRING database to define key protein-protein interactions (Yu et al., 2012). The network was constructed using Cytoscape 3.6.1 (Shannon et al., 2003).

Validation of RNA-seq by quantitative RT-PCR (qRT-PCR)
To confirm the RNA-seq data, 20 DEGs (10 with increased expression and 10 with decreased expression) were selected for qRT-PCR validations. Reactions were performed using SYBR Premix Ex Taq TM Kit (Takara) according to the manufacturer's protocol. Reaction conditions were as follows: 40 cycles of 95 • C for 30 s; 95 • C for 5 s; and 60 • C for 30 s. PCRs were performed on a BIO-RAD CFX-Connect Real-Time System. Relative gene expression was determined using the 2 − Ct method normalized to 18S rRNA (Livak & Schmittgen, 2001). Significant differences were determined using a t -test with a threshold of p < 0.05. Primers are listed in Table S1. Primer efficiency and melting curves are listed in Table S2 and Figs. S1-S5.

MIC assays and time-kill curves
The MIC of MAF-1A against C. parapsilosis was determined as 0.6 mg/mL. Time-kill curves of MAF-1A at MIC showed a gradual antifungal effect during the first 8 h of C. parapsilosis culture (Fig. 1). After 8 h, cell numbers increased but remained lower than those of the control group.

Transcriptional stress responses and enrichment analysis of MAF-1A treated C. parapsilosis
RNA-seq analysis in C. parapsilosis treated with MAF-1A for 6 and 18 h showed 5,747 DEGs. Sequence reads were deposited in the NCBI Sequence Read Archive (SRA) under accession number PRJNA638006. A total of 2,439 DEGs were detected. Out of these genes, 2187 were identified at 6 h, representing 38.05% of the total detectable genes. A total of 252 genes were differentially expressed after 18 h and accounted for 4.38% of the total expressed genes. After 6 and 18 h of MAF-1A treatment, 67 DEGs with opposite trends were observed (reversed genes 1, RG1 and reversed genes 2, RG2). In total, 56 DEGs were upregulated, whilst down-regulated genes remained unchanged (one unchanged genes: UG1; two unchanged genes: UG2) (Fig. 2).

DEG enrichment analysis
Amongst the DEGs at 6 h, 1122 showing increased expression were enriched in 85 KEGG pathways, 20 of which were significant with padj < 0.05. The most significant pathways with increased expression in C. parapsilosis following MAF-1A treatment were: oxidative phosphorylation, peroxisome, citrate cycle (TCA cycle), carbon metabolism, cell cycleyeast, MAPK signaling-yeast, meiosis-yeast, fatty acid degradation, and autophagy ( Table 1). Genes of decreased expression were enriched in steroid biosynthesis, biosynthesis of amino acids, cysteine and methionine metabolism, biosynthesis of antibiotics, ribosome, RNA polymerase, biosynthesis of secondary metabolites, RNA transport, ribosome biogenesis in eukaryotes, lysine biosynthesis, 2-oxocarboxylic acid metabolism, pyrimidine metabolism, glycine, serine and threonine metabolism and purine metabolism (Table 2). At 18 h, 101 genes were upregulated and enriched in 24 KEGG pathways, two of which were significant (Table 3). A total of 151 genes were downregulated and significantly enriched in carbon metabolism, biosynthesis of antibiotics, oxidative phosphorylation, the biosynthesis of secondary metabolites, and the biosynthesis of amino acids (Table 4).

RG1, RG2, UG1, and UG2 gene enrichment analysis
The 42 genes in RG1 were enriched in 17 KEGG pathways, of which oxidative phosphorylation was most significant. Additionally, 25 genes in RG2 were enriched in 13 KEGG pathways, of which arginine biosynthesis, the biosynthesis of antibiotics, the biosynthesis of amino acids, and the biosynthesis of secondary metabolites were enriched. A total of 24 genes in UG1 and 32 genes in UG2 were enriched in 8 and 20 KEGG pathways, respectively. In UG1, the genes were enriched in butanoate metabolism, propionate metabolism, beta-alanine metabolism, valine, leucine, and isoleucine degradation. No pathways were significantly enriched in UG2 at a padj < 0.05 (Table 5). Genes in RG1, RG2, UG1, and UG2 were enriched in 535 GO terms, a total of 9 of which were significant (Table 6). Genes in RG1 were involved in energy production and redox processes. Genes in RG2 were associated with the anabolic processes of various organic acids. Genes in UG2 were involved in oxidation-reduction processes.

Verification of differentially expressed genes
A total of 20 genes were selected, including 10 with increased expression and 10 with decreased expression. Genes were evenly selected from 6 h and 18 h time points to validate the RNA-seq data by qRT-PCR. The results indicated the expression levels have a  consistent change for both RNASeq and qRT-PCR. Hence, the qRT-PCR results confirmed the reliability of our RNA-Seq data (Fig. 3).

DISCUSSION
C. parapsilosis is one of the most prevalent fungal species in many regions. In addition to its high rates of infection, its etiology differs from that of C. albicans (Holland et al., 2014). Specific C. parapsilosis isolates are resistant to conventional antifungal drugs including echinocandins, azoles, and amphotericin B (Lotfali et al., 2016;Maria et al., 2018;Thomaz et al., 2018). Antimicrobial peptides lead to cell lysis and death through  cell membrane leakage (Papo & Shai, 2003;Paterson et al., 2017;Shai, 1999;Utesch et al., 2018). However, mechanistic studies of antimicrobial peptides have determined that their membrane interactions are complex. Park, Kim & Kim (1998) showed that buforin II prevents microorganisms entry into cells. Lee et al. (2019) found that antifungal β-peptides cause cell death by entering cells and causing nuclear and vacuole dysfunction. Chileveru et al. (2015) showed that human alpha-defensin 5 enters the cytoplasm of Escherichia coli and interferes with cell division. AMPs work through various mechanisms. In our previous studies, we showed that MAF-1A inhibits C. albicans through its effects on the cell wall, cell membrane, and ribosomes (Wang et al., 2017). In this study, we found that MAF-1A alters gene expression in several important biological pathways in C. parapsilosis, including oxidation-reduction processes and alternative energy sources. We further compared the response of C. albicans and C. parapsilosis to MAF-1A, most DEGs have the same expression trend (upregulated/downregulated), and identified changes in both stress and energy metabolism pathways (carbon metabolism, cell cycle, peroxisome, carbon metabolism, fatty acid degradation) (Tables S3 and S4). We hypothesized that the antifungal peptide MAF-1A exerts antifungal effect and disrupts energy metabolism by affecting oxidation-reduction processes, due to its effects on the mitochondria. Whilst antimicrobial peptides have multiple modes of action, these remain undetermined for MAF-1A. Our findings suggest that intracellular targets may be the key sites of MAF-1A activity, with enrichment analysis of the DEGs suggesting that MAF-1A exerts antimicrobial activity through a variety of mechanisms.

Mitochondrial function
A total of 42 genes showed increased expression at 6 h and were enriched in oxidative phosphorylation pathways. The genes were mainly involved in respiratory chain electron transport processes in the mitochondrial inner membrane, including NHD release of H+ and ATP synthase (Fig. 7) (Kanehisa et al., 2019;Ogata et al., 1999). A total of six genes in RG1 were also enriched in this pathway (CpCOX15,CpQCR2,CpQCR8,CpQCR7,and CpCOR1). Previous studies showed that COX15 encodes an indispensable mitochondrial protein for Saccharomyces cerevisiae cytochrome oxidase (Glerum et al., 1997). Cytochrome oxidase is a terminal enzyme in the respiratory electron transport chain that is essential for ATP synthesis. Reactive oxygen species (ROS) are produced by the oxidative phosphorylation of ATP and can disrupt the electron transport chain in mitochondria (Piippo et al., 2018). ROS production induces damage to lipids, proteins, lipids, and nucleic acids, leading to cell death. Eukaryotes prevent cell damage through oxidative stress detoxification and the prevention of ROS accumulation. In this study, GO enrichment analysis of the RG1 genes showed that 7 that were highly expressed were associated with redox processes (CpPOX9, CpCOX15, Cpnuo-21, CpAAEL001134, CpQCR7, CpRGI1, and CpnamA). These processes help cells to remove accumulated ROS. Due to the increased expression of these genes in response to MAF-1A, it is possible that MAF-1A promotes oxidative phosphorylation, which disrupts electron transfer in the mitochondria, enhancing ROS production and subsequent cell damage. C.parapsilosis has an unusual mitochondrial genome architecture, consisting of linear DNA molecules of 30.9-Kbp, terminating with specific telomeric structures on both sides (738-Kbp long). This differs from telomeres at the ends of eukaryotic nuclear chromosomes, particularly in humans (Kovac, Lazowska & Slonimski, 1984). MAF-1A interferes with the expression of multiple genes related to the mitochondrial functions of C. parapsilosis, It is, therefore, feasible that MAF-1A interferes with the normal function of C. parapsilosis mitochondria.

CONCLUSIONS
In summary, MAF-1A has a complex response in C. parapsilosis. Most DEGs identified through RNA-seq analysis were related to oxidation-reduction processes and the use of alternative energy sources, Mitochondria are important target for the anti-fungal peptide MAF-1A to exert anti-C. parapsilosis. RNA-seq data therefore provide future direction to study the antifungal mechanisms of MAF-1A and highlight the potential pathways that contribute to resistance.

Grant Disclosures
The following grant information was disclosed by the authors: