A Novel Butanol Tolerance-Promoting Function of the Transcription Factor Rob in Escherichia coli

Producing high concentrations of biobutanol is challenging, primarily because of the toxicity of butanol toward cells. In our previous study, several butanol tolerance-promoting genes were identified from butanol-tolerant Escherichia coli mutants and inactivation of the transcriptional regulator factor Rob was shown to improve butanol tolerance. Here, the butanol tolerance characteristics and mechanism regulated by inactivated Rob are investigated. Comparative transcriptome analysis of strain DTrob, with a truncated rob in the genome, and the control BW25113 revealed 285 differentially expressed genes (DEGs) to be associated with butanol tolerance and categorized as having transport, localization, and oxidoreductase activities. Expression of 25 DEGs representing different functional categories was analyzed by quantitative reverse transcription PCR (qRT-PCR) to assess the reliability of the RNA-Seq data, and 92% of the genes showed the same expression trend. Based on functional complementation experiments of key DEGs, deletions of glgS and yibT increased the butanol tolerance of E. coli, whereas overexpression of fadB resulted in increased cell density and a slight increase in butanol tolerance. A metabolic network analysis of these DEGs revealed that six genes (fadA, fadB, fadD, fadL, poxB, and acs) associated with acetyl-CoA production were significantly upregulated in DTrob, suggesting that Rob inactivation might enhance butanol tolerance by increasing acetyl-CoA. Interestingly, DTrob produced more acetate in response to butanol stress than the wild-type strain, resulting in the upregulation expression of some genes involved in acetate metabolism. Altogether, the results of this study reveal the mechanism underlying increased butanol tolerance in E. coli regulated by Rob inactivation.

Introduction

Butanol has received increased attention as a high-energy-density fuel and bulk chemical feedstock (Gu et al., 2011). Only biobutanol can be used in the pharmaceutical and cosmetics industries, and biomass-based butanol production can reduce environmental pollution compared to petrochemical-based production (Jiang et al., 2015b). Traditional biobutanol is produced by Clostridium species and can reach a titer of approximately 20 g/L (Qureshi and Blaschek, 2000). The butanol titers produced by engineered Escherichia coli strains can also reach approximately 20 g/L in batch fermentation without any antibiotics or inducers and up to 30 g/L with the continuous removal of butanol (Shen et al., 2011; Jang and Lee, 2015; Dong et al., 2017), indicating that E. coli is a potential butanol-producing bacterium. However, the current production of butanol is limited to approximately 20 g/L due to the toxicity of butanol toward bacteria; indeed, this toxicity is a bottleneck in butanol production (Qureshi and Blaschek, 2000). Therefore, it is important to improve the butanol tolerance of chassis strains and to explore tolerance-associated mechanisms to promote the highly efficient production of butanol.

Due to its hydrophobicity, butanol binds to lipid chains of the bacterial cytoplasmic membrane, which results in toxicity to the cells (Jones and Woods, 1986; Nielsen and Prather, 2009). The mechanism of this toxicity involves increased membrane fluidity and permeability in the presence of butanol, causing leakage of protons and ATP and interference with the correct folding of proteins, which leads to cell damage or even death (Sikkema et al., 1995; Isken and de Bont, 1998). Bacteria alter their physiological and biochemical characteristics to respond to butanol stress in the following ways: (1) membrane composition changes, including increasing efflux capacity of toxins, preventing leakage of intracellular components, and protecting cells from damage due to the solvent (Reyes et al., 2011; Bui et al., 2015; Royce et al., 2015; Sandoval and Papoutsakis, 2016); (2) physiological responses similar to those involved in responses to osmotic, oxidizing, respiratory, and heat shock stresses (e.g., altered osmotic pressure) (Purvis et al., 2005; Chin et al., 2017), reactive oxygen species (ROS) accumulation (Kaczmarzyk et al., 2014) and enhanced metabolic transport and molecular chaperone levels (Peralta-Yahya et al., 2012); and (3) up- or downregulation of expression of regulatory genes, such as those encoding sensor proteins, transcription factors (TFs), and those involved in regulating expression of small RNAs (Peralta-Yahya et al., 2012), to modulate gene expression profiles to protect against butanol stress (Reyes et al., 2011, 2012).

Escherichia coli expresses 304 TFs (Perez-Rueda et al., 2015), though only a few have been reported to regulate butanol tolerance-related genes (Reyes et al., 2012; Horinouchi et al., 2018). Furthermore, the corresponding genes regulated by these TFs and the associated regulatory mechanism have not yet been clarified (Aquino et al., 2017), limiting the improvement of butanol-tolerant chassis strains using a rational-design engineering strategy. In our previous study, a mutant strain (BW1847) able to tolerate 2% (v/v) butanol was obtained, and among the genes mutated or deleted, rob (GenBank No. RS22900), acrB (GenBank No. RS02385), and tqsA (GenBank No. RS08380) have been identified as having a function in enhancing butanol tolerance (He et al., 2019). The DTrob strain, with an AT_686–7 base deletion within rob in the genome, and the rob-deletion mutant Δrob produced much higher cell densities than did the wild-type strain under 0.75% butanol stress (He et al., 2019), indicating that partial and full inactivation of Rob both result in improved butanol tolerance. The rob gene encodes a right oriC-binding transcriptional activator (Nakajima et al., 1995) that interacts with a superoxide response regulon transcriptional activator (SoxS) and a multiple antibiotic resistance transcriptional regulator (MarA) and can activate genes involved in antibiotic, oxygen pressure and organic solvent resistance (Nakajima et al., 1995). The rob-deleted strain exhibited decreased tolerance to antibiotics, oxygen stress, cyclohexane and n-pentane (Nakajima et al., 1995; White et al., 1997; Bennik et al., 2000). However, inactivation of Rob yielded increased butanol tolerance in our previous study. These contrasting results are due to the extensive, intricate and multiple regulatory mechanisms of Rob, and the corresponding butanol-tolerant mechanism regulated by Rob is unknown. Strain DTrob with truncated Rob exhibited a butanol tolerance characteristic and can be used to investigate the interaction mechanism between Rob and its target genes. DTrob is thus an ideal candidate for studying the butanol tolerance mechanism caused by Rob inactivation.

Therefore, in this study, butanol stress response genes regulated by the Rob-inactivated mutant were evaluated by RNA-Seq. In addition, key genes involved in butanol tolerance were functionally identified, and a potential tolerance mechanism was determined to demonstrate the novel roles of Rob in response to butanol stress.

Materials and Methods

Growth Assays and Extracellular n-Butanol Measurements

The site-specific mutant DTrob (AT_686–7 deletion in the rob gene) was constructed in a previous study (He et al., 2019). The DTrob strain and BW25113 (Supplementary Table S1) were precultured in lysogeny broth (LB) until the late exponential phase, after which the cultures were concentrated to an OD₆₀₀ of 20. The concentrated cells were then used to inoculate 50 mL of LB medium containing 0, 0.75, 1, and 1.25% (v/v) butanol at an initial OD₆₀₀ of 0.1–0.15 to evaluate cell growth. A total of 0.7 mL of culture was taken to measure the OD₆₀₀ with a Cary 50 Conc spectrophotometer (Varian, Palo Alto, CA, United States); 1 mL of culture was collected at the appropriate time to measure butanol and acetate concentrations in the medium. The collected samples were centrifuged, and the supernatants were used to measure butanol and acetate contents by high-performance liquid chromatography (HPLC). HPLC was performed using an organic acid analytical column (Aminex HPX-87H Ion Exclusion Column, 300 mm × 7.8 mm) at 45°C; sulfuric acid (4 mM) was used as the mobile phase at 0.8 mL/min. The butanol and acetate concentrations in each sample were calculated by comparisons with the peak area of the standard. The per unit intracellular butanol concentration (PIC) was calculated according to previous reports (He et al., 2019), with a slight modification, as shown by the following equation:

$$\text{PIC}=\frac{{\text{C}}_{\text{initial}}-{\text{C}}_{\text{final}}}{\text{N}}$$

where “PIC” indicates “per unit intracellular butanol concentration” (μg/L) and “C_initial” and “C_final” indicate “initial extracellular butanol concentration” and “final extracellular butanol concentration”, respectively; “N” indicates the number of cells; and 1OD cells at 600 nm corresponds to 8.3 × 10⁸ cells mL^–1. N = OD₆₀₀ × 8.3 × 10⁸ mL^–1. The values of the mean and standard deviation are plotted using the bar and error bar.

Overexpression strains (Supplementary Table S1) were precultured in LB until the late exponential phase, and the cultures were concentrated to an OD₆₀₀ of 20. The concentrated cells were then used to inoculate LB medium containing 0.2% (w/v) L-arabinose, 100 μg/mL ampicillin, and a gradient concentration of butanol at an initial OD₆₀₀ = 0.2. A total of 1 mL of culture was collected at a suitable time point for cell concentration measurements.

RNA-Seq and Data Processing

The strains BW25113 and DTrob were cultured in 50 mL of LB containing 0.75% (v/v) butanol until the OD₆₀₀ reached 0.6–1.1; 3-mL aliquots were harvested, and the cells were washed with ice-cold PBS (phosphate buffer saline) buffer. Total RNA was extracted with TRIzol reagent (Invitrogen, Carlsbad, CA, United States) according to the manufacturer’s instructions. The quality of the total RNA extracted was assessed spectrophotometrically at 230, 260, and 280 nm. The RNA sample with an OD₂₆₀/OD₂₃₀ ratio higher than 1.8 and an OD₂₆₀/OD₂₈₀ ratio between 1.8 and 2.1 was considered as pure and used to cDNA library construction. cDNA library preparation and sequencing were conducted by the Allwegene Technology Company in Beijing, China. Individual libraries were sequenced using the Illumina HiSeq 4000 platform (Illumina, Inc., San Diego, CA, United States). The original fluorescent images were converted into raw sequence reads by CASAVA software. The raw reads were then processed to remove reads with adapters, reads containing poly-N (N > 1%) and low-quality reads (more than 50% bases with Q-score ≤20), yielding clean reads (Gao et al., 2018; Wang et al., 2018). The clean data were deposited in Gene Expression Omnibus (GEO) at NCBI (Accession Number GSE120032) and mapped to the genome of the reference strain BW25113 (LOCUS NZ_CP009273) using Bowtie2 (Kim et al., 2018). Gene expression levels were analyzed using fragments per kilobase of exon model per million mapped reads (FPKM). Differentially expressed genes (DEGs) between the two samples were identified using the DESeq package (v1.24.0), and | log2(FoldChange)| > 1 and q-value < 0.005 were used as the threshold to identify significant differences in gene expression. Analysis of Gene Ontology (GO), which assigns genes into functional categories, was performed using the GOSeq (v1.22) and top GO (v2.22) R packages.

Quantitative Reverse Transcription PCR Analysis

The expression levels of eighteen DEGs representing different functional categories were assessed by real-time quantitative reverse transcription PCR (qRT-PCR) to validate the reliability of the RNA-Seq data. Total RNA from DTrob and BW25113 cells was extracted as described above, and RQ1 RNase-Free DNase (Promega, Madison, WI, United States) was added to the RNA to remove genomic DNA. PCR reaction was performed using the total RNA as template in order to check for genomic DNA contamination. No product was observed, demonstrating that the sample had no genomic DNA contaminants and was suitable for qRT-PCR assay. Subsequently, reverse transcription was performed using iScript cDNA Synthesis Kit (Bio-Rad) with 2 μg of total RNA following the manufacturer’s instructions. qRT-PCR experiments were performed using CFX96 Real-Time System (Bio-Rad). Each reaction contained 2 μL of diluted (1/10) cDNA, Taq SYBR Green qPCR Premix (Yugong Biolabs Inc., Jiangsu, China) and the corresponding primer pairs (Supplementary Table S2). The 16S rRNAgene was used as a reference to normalize the qRT-PCR data. The expression level of each gene was calculated according to the following formula: expression level = 2^–ΔΔCt (Livak and Schmittgen, 2001). The relative mRNA level is presented as the percent (%) ratio of the gene expression level in DTrob to that in BW25113.

Construction of Knockout Strains

Strains deleted for the glgS and yibT genes were constructed using the CRISPR-Cas9 system (Ran et al., 2013; Jiang et al., 2015a; Supplementary Figure S1). Inverse PCR was performed using pTargetF plasmid as a template to introduce the target sequence of N20 (20-bp complementary region) upstream of the single-guide RNA (sgRNA) in the pTargetF plasmid. The N20 sequence was introduced into the primers shown in Supplementary Table S3. The inverse PCR products were digested with the methylation-sensitive restriction enzyme DMT (TransGen Biotech, Shanghai, China) to remove methylated plasmid templates, after which the DNA was transformed into DH5α competent cells (Biomed, Beijing, China). Positive clones were identified by PCR with the appropriate geneN20F/pTargetF-IR primer pairs (Supplementary Table S3), and the resulting plasmid was correspondingly named pTargetF-geneN20. DNA fragments containing left and right homologous arms located upstream and downstream of the target gene were amplified with the primer pairs geneDLF/geneDLR and geneDRF/geneDRR (Supplementary Table S3), respectively. The two fragments were then fused together by overlap PCR with the appropriate geneDLF/geneDRR primer pairs (Supplementary Table S3), and the resulting fused PCR fragment was used as donor DNA to delete the target gene by homologous recombination. Approximately 500 ng of donor DNA and the corresponding pTargetF-geneN20 plasmid were co-transformed into BW25113 (pCas) competent cells. Positive clones were screened on LB plates containing 100 mg/L kanamycin and spectinomycin and subsequently identified by PCR using the appropriate geneDLF/geneDRR primer pairs (Supplementary Table S3). The deletion was confirmed by sequencing, and the positive clones were cultured in LB medium supplemented with 50 mg/L kanamycin and 0.5 mM IPTG to eliminate the pTargetF-geneN20 plasmid. The temperature-sensitive plasmid pCas was then removed by growing the culture overnight at 37°C.

Overexpression Analysis of Several Upregulated Genes

Thirteen upregulated genes (Supplementary Table S5) were chosen for overexpression analysis. Plasmids harboring the target genes were constructed on the basis of pBAD30 (Guzman et al., 1995). PCR amplification of the target genes from the genomic DNA of BW25113 was performed using the primer pairs shown in Supplementary Table S4. The plasmid pBAD30 was digested with the restriction endonucleases SacI and HindIII, and the aforementioned PCR products were cloned into the SacI/HindIII sites of pBAD30 using GenBuilder Cloning Kit (GenScript, Nanjing, China). The mixture was transformed into competent DH5α cells and screened on LB agar plates containing 100 μg/mL ampicillin. Single clones were identified by PCR using the primers pBADIup/pBADIdown (Supplementary Table S4). Plasmids were extracted using TIANprep Mini Plasmid Kit (TIANGEN, Beijing, China) and verified by sequencing (Sangon Biotech, Shanghai, China). The overexpression plasmids were transformed into BW25113 by electroporation, and the positive transformants were correspondingly named BW25113 (gene) (Supplementary Table S1).

Results

Enhanced Butanol Tolerance by the Truncated Rob in DTrob

The DTrob strain has an AT_686–7 base deletion of rob in the genome, which leads to an early termination of translation of rob mRNA, and the resulting truncated protein only has 229 amino acids (He et al., 2019). Strain DTrob showed a similar growth trend to control BW25113 in absence of butanol (Figure 1A), and the maximum cell density of the DTrob exhibit 94, 49, and 18% higher than that of BW25113 under 0.75, 1, and 1.25% (v/v) butanol stress, respectively, indicating that the DTrob strain was able to tolerate 1–1.25% and that the improvement in relative growth decreased with increasing butanol concentrations. The PIC of the DTrob strain was 76, 56, and 78% of that of the control BW25113 under 0.75, 1, and 1.25% (v/v) butanol stress at 4 h (logarithmic phase), respectively (Figures 1B–D), which shows that inactivation of Rob can result in an enhanced ability to efflux butanol out of the cell to improve the tolerance, and the highest efflux capacity was shown in the presence of 1% butanol. The PIC value of DTrob was 32–60% of that of BW25113 in 0.75–1% (v/v) butanol at 12 h, demonstrating a higher efflux capacity than at 4 h. Nevertheless, interestingly, the PICs of both BW25113 and DTrob in a higher concentration of butanol (1.25%) were decreased at 12 h compared to 4 h, and the two strains had similar intracellular butanol concentration. It is due to that the strains grown in 1.25% butanol have entered the decline phase at 12 h (Figure 1D), and the cell death causes more butanol to be released into the medium, thus resulting in a decreased sharply butanol content in these cells and lower PIC value than at 4 h. Therefore, the alleviated toxicity by butanol efflux in log stage is one reason for the improved butanol tolerance of DTrob.

FIGURE 1

Figure 1. Growth assays for the DTrob strain in the presence of butanol. (A) Without butanol; (B) with 0.75% (v/v) butanol; (C) with 1% butanol; and (D) with 1.25% butanol. The line graph shows the growth curve, and the bar graph shows the unit intracellular butanol concentration. Error bars represent one standard deviation. Statistical significance of the differences between DTrob and BW25113 was analyzed by a one-way analysis of variance (ANOVA). *P < 0.05 was considered statistically significant and **p < 0.01 was considered extremely significant. Growth assays were performed in triplicate in screw-cap flasks.

Transcriptomic Differences Between the DTrob and Control Strains Under Butanol Stress

RNA-Seq of the wild-type BW25113 and DTrob strains was performed to assess differences in gene expression between them. Genes with a “q-value < 0.005” and “| log2(FoldChange)| > 1” were defined as DEGs. A total of 285 (6.2% of genes) DEGs were identified in DTrob compared with the control. In the late logarithmic phase, 184 and 101 genes were upregulated and downregulated in DTrob compared with BW25113, respectively (Supplementary Figure S2A). The 285 DEGs were then subjected to GO enrichment analysis and grouped into three primary categories: biological process, molecular function and cellular component. GO annotations were found for 202 of these genes. The significantly enriched GO terms and number of DEGs with GO annotations were as follows: “localization” (32.7%, 66/202), “transport” (31.7%, 64/202), “establishment of localization” (31.7%, 64/202), “oxidoreductase activity” (17.82%, 36/202), and “transporter activity” (20.30%, 41/202) (Supplementary Figure S2B). In addition, the 285 DEGs were examined by Kyoto Encyclopedia of Genes and Genomes (KEGG) pathway analysis using KOBAS (v2.0), revealing 68 pathways, 53 of which were different from the KEGG pathways in previous reports (Reyes et al., 2012; Si et al., 2016; Guo et al., 2019). The significant KEGG terms are primarily associated with the biosynthesis of antibiotics, biosynthesis of secondary metabolites and fatty acid degradation (Supplementary Figure S3).

Hundreds of potential genes involved in the n-butanol tolerance response have been identified by DNA microarrays and comparative genome hybridization microarrays (array-CGH); nevertheless, only approximately 20 candidate genes were identified as being associated with butanol tolerance through overexpression or knockout experiments (Supplementary Table S6; Rutherford et al., 2010; Reyes et al., 2011; Si et al., 2016). Among the 285 DEGs identified by RNA-Seq in the present study, ompF, acrB, glcF, ybjC, yibT and cpxP have been reported to be associated with butanol tolerance (Rutherford et al., 2010; Reyes et al., 2011; Si et al., 2016); however, most of the DEGs have not yet been reported to be involved in butanol tolerance. Therefore, the results of this study reveal potential functional genes related to the development of butanol-tolerant strains for further exploration of butanol tolerance mechanisms.

Validation of RNA-Seq Data Using Real-Time qRT-PCR

Twenty five genes representing different GO functional categories and different expression levels were classified based on annotations in the NCBI database (Table 1) and primarily classified into four groups: regulatory factors (ChaC, GadE, GlgS, and YiaG), transport and membrane proteins (ActP, ElaB, PspG, FadL, MdtD, OmpW, Slp, YacH, YibT, YdcL, and YjcH), stress response proteins (HdeB, InaA, and YhbO), and enzymes (Acs, AstC, NarG, FadA, FadB, FadD, and PoxB). The expression levels of 25 DEGs between the DTrob and BW25113 strains under butanol stress were analyzed by qRT-PCR to confirm the reliability of the RNA-Seq data obtained. The log2-fold change values of the 25 DEGs were also examined to compare the consistency of the RNA-Seq and qRT-PCR data, and a strong correlation (R = 0.8318) was observed between the RNA-Seq and qRT-PCR data (Supplementary Figure S4). Based on the qRT-PCR results, the expression trends of 23 of these candidate genes were consistent with the RNA-Seq data. Among these 23 genes, acs, actP, astC, elaB, fadA, fadB, fadD, fadL, yiaG, yhbO, and yjcH were obviously upregulated (263–11503%) (Table 1) and gadE, hdeB, and slp slightly upregulated (153–194%) in the DTrob strain compared to the control strain, whereas chaC, glgS, inaA, narG, mdtD, yibT, ydcL, and yacH were notably downregulated (65–99%) (Table 1). Therefore, 92% of the candidate genes exhibited consistent qRT-PCR and RNA-Seq results, indicating that the RNA-Seq data were of good quality.

TABLE 1

Table 1. Functional categories and qRT-PCR results of 25 DEGs associated with butanol tolerance.

Functional Analysis of Significant Differentially Expressed Genes

Based on GO enrichment and the gene expression levels of the DEGs, fifteen representative DEGs in each cluster were selected (Supplementary Table S5) for investigation of their function in promoting butanol tolerance by growth evaluation under butanol stress. To determine their direct impacts on the butanol tolerance phenotype, the selected upregulated and downregulated genes were overexpressed and deleted, respectively. Deletion mutants of glgS and yibT were obtained using the CRISPR-Cas9 system and named ΔglgS and ΔyibT, respectively. The growth evaluation results showed that the maximum cell density of the ΔglgS (pBAD30) and ΔyibT (pBAD30) strains was increased by 7.7–23.6% compared to that of BW25113 in the presence of 1–1.25% (v/v) n-butanol, indicating that the absence of Rob or downregulated expression improved butanol tolerance (Figure 2). The complementation strains ΔglgS (pBAD-glgS) and ΔyibT (pBAD-yibT) of the two genes displayed a lower cell density than their deletion strains ΔglgS and ΔyibT, indicating that their overexpression can lead to a decrease in butanol tolerance (Figures 2C,D). Overexpression of these two genes in the wild-type strain also reduced butanol tolerance compared to the control strain BW25113 (pBAD30) (Figure 2D). These results indicate that a reduced level of glgS and yibT expression is beneficial for enhancing butanol tolerance.

FIGURE 2

Figure 2. Butanol tolerance evaluation of the glgS or yibT deletion and complement strains. All strains were cultured in LB medium containing 0% (A), 0.75% (B), 1% (C), and 1.25% (v/v) (D) butanol.

Thirteen upregulated genes verified by qRT-PCR were selected for investigating their function in promoting butanol tolerance, and engineered strains were cultured in the presence or absence of butanol. The maximum cell density of BW25113 (pBAD-fadB) was increased by 25.5, 12.7, and 5.1% under 0, 0.5, and 0.75% (v/v) butanol stress, respectively, compared with that of the control strain BW25113 (pBAD30) (Figure 3). These results indicate that fadB overexpression can promote cell growth and that the cell density of E. coli under butanol stress can also be increased to some extent. The maximum cell density of the slp-, acs-, and fadL-overexpressing strains was decreased by 31.9, 8.3, and 27.7%, respectively, compared with the control strain in the absence of butanol (Figure 3A, with a significant decrease under butanol stress (Figures 3B,C). These results indicate that overexpression of these three genes inhibits the growth of the cells, resulting in a decreased cell density under butanol stress.

FIGURE 3

Figure 3. Growth evaluation of overexpression strains. (A–C) Growth curves of acs-, fadB-, fadL-, and slp-overexpressing strains exposed to 0, 0.5, and 0.75% (v/v) butanol. BW25113 (pBAD30) was used as a control.

Unlike fadB, slp, acs, and fadL, overexpression of the other nine upregulated genes did not cause an increase or decrease in cell growth under butanol stress (data not shown), indicating that overexpression of these nine genes does not alter butanol tolerance. It is possible that most responsive genes do not exhibit a significant butanol tolerance function; it also may be that simultaneous changes in the expression of multiple DEGs cause synergistic physiological and biochemical responses to resist butanol stress. Therefore, the function and classification of the DEGs is not well illustrated merely by functional complementation via knockout or overexpression of a certain gene. Correspondingly, the metabolic network of these DEGs was also analyzed, and the following interesting response mechanism was observed.

DTrob Responds to Butanol Stress by Altering Acetyl-CoA and Acetate Production

The expression levels of poxB, encoding pyruvate dehydrogenase (PoxB), and acs, encoding acetyl-CoA synthetase (Acs), were upregulated in the DTrob strain, as shown by the RNA-Seq and qRT-PCR data (Figure 5). PoxB and Acs convert pyruvate to acetyl-CoA in a two-step enzymatic reaction (Figure 4A). Four genes (fadA, fadB, fadD, and fadL) involved in fatty acid β-oxidation were also significantly upregulated in the DTrob strain, and the final product of fatty acid β-oxidation is acetyl-CoA (Figure 4A). An enhanced acetyl-CoA pool has been demonstrated to be achieved through modification of a fatty acid β-oxidation pathway (Zhang et al., 2019), and acetyl-CoA is a key molecule in microbial central carbon metabolism and involved in a variety of cellular processes (Krivoruchko et al., 2015). Thus, upregulation of the six genes (Figure 5 and Table 1) caused by inactivation of Rob might result in increased acetyl-CoA levels, which would provide more raw material and energy for a number of physiological processes and cell growth. Nonetheless, the acetyl-CoA content of BW25113 and DTrob was not significantly different under 0.75% butanol stress (data not shown). It is possible that the increased acetyl-CoA in DTrob is quickly transferred to other metabolic pathways, such as the TCA cycle to provide energy for cell growth, or provides acetyl groups for acetylation of histones to promote cell growth and proliferation (Cai et al., 2011). Accordingly, no obvious accumulation of acetyl-CoA was observed in DTrob, although which showed a higher cell density than BW25113 under butanol stress.

FIGURE 4

Figure 4. Mechanism of butanol tolerance regulated by Rob via an increase in acetyl-CoA. (A) The fatty acid oxidation and acetate metabolism pathways regulated by Rob. Genes marked in red font were upregulated in DTrob. (B) Acetate levels in the BW25113 and DTrob strains under 0.75% (v/v) butanol stress. The line graph shows the growth curve, and the bar graph shows the extracellular acetate concentration. Error bars represent one standard deviation. Statistical significance of the differences between DTrob and BW25113 was analyzed by a one-way ANOVA. *P < 0.05 was considered statistically significant and **p < 0.01 was considered extremely significant. Growth assays were carried out in triplicate in screw-cap flasks.

FIGURE 5

Figure 5. qRT-PCR analysis of genes involved in acetyl-CoA production. Gray bar represents BW25113, blue bar represents DTrob. The expression level was calculated by 2^–ΔΔCt method.

In addition, expression of ten genes (Slp, hdeA, hdeB, hdeD, gadA, gadB, gadC, gadE, yiaG, and yccJ) upregulated in the DTrob strain has been reported to be increased when cells are subjected to acid stress (Ma et al., 2003; Tucker et al., 2003; Hommais et al., 2004). Indeed, upregulation of PoxB, which converts pyruvate to acetate (Dittrich et al., 2005), would lead to an increase in acetate levels, and the extracellular acetate concentration in the DTrob culture was higher than that of the BW25113 strain in the log phase (Figure 4B). These results indicate that upregulated expression of some genes in the DTrob strain might promote the production of acetyl-CoA and result in acetate accumulation, in turn inducing expression of genes involved in the response to acid stress.

Discussion

Mutation of the transcription regulator rob in E. coli altered the expression levels of 285 genes in response to butanol stress. Deletion of the genes glgS and yibT, encoding a transcriptional regulator and a membrane-associated protein, respectively, directly increased butanol tolerance. In addition, overexpression of the fadB gene, encoding a fatty acid oxidation-related enzyme, promoted the growth of E. coli and slightly enhanced butanol tolerance. Deletion of glgS has been reported to improve the production of both flagella and type 1 fimbriae, resulting in enhanced swarming motility of bacterial cells. Such enhanced motility protects cells from toxicity by allowing them to escape from toxic environments, such as those containing n-butanol (Rahimpour et al., 2013; Kim et al., 2016). However, the ΔglgS strain showed an increased capacity to initiate biofilm formation (Rahimpour et al., 2013), and biofilms provide good protection for cells as a permeability barrier and also enhance the ability of bacteria to adapt to adverse environments, including butanol stress (Li et al., 2016). It has been reported that deletion of yibT can increase the proportion of unsaturated fatty acids in the membrane, resulting in lower membrane fluidity and higher membrane rigidity and integrity, protecting cells from organic solvent stress (Si et al., 2016). Therefore, yibT deletion can improve the butanol tolerance of E. coli.

The gene fadB encodes the α-subunit of a multienzyme complex that is involved in fatty acid β-oxidation (Fujita et al., 2007) and responsible for the hydration of enoyl-CoA and oxidation of 3-hydroxyacyl-CoA (Figure 4A). Overexpression of fadB may accelerate the two-step enzymatic reaction and the conversion of fatty acids to acetyl-CoA. Additional acetyl-CoA entering the tricarboxylic acid cycle would provide more energy for cell growth, resulting in an increase in cell density. This growth improvement may also result in a higher cell density under butanol stress.

Interestingly, overexpression of acs, slp, and fadL did not increase the butanol tolerance of E. coli and even led to inhibition of cell growth. The acs, slp, and fadL genes encode the Acs (Brown et al., 1977), outer membrane protein Slp (Alexander and St John, 1994) and fatty acid transporter FadL (van den Berg et al., 2004), respectively, and the results of previous studies also showed that overexpression of these three genes inhibits cell growth (Black et al., 1985; Alexander and St John, 1994; Zhang et al., 2019). However, individually overexpressing these three genes does not improve butanol tolerance, as they function in synergy with other genes in the same GO functional category to resist butanol stress. Moreover, Acs and PoxB may function together to promote the production of acetyl-CoA; FadL may act in conjunction with proteins encoded by fadA, fadB, and fadD to promote fatty acid degradation, and Slp may act in conjunction with some proteins to promote acid tolerance in response to acetate stress (Figure 4A). Therefore, the proteins encoded by these genes synergistically respond to butanol stress.

Among the 15 DEGs identified by functional complementation experiments, only fadB, glgS, and yibT were identified to be involved in significant butanol tolerance in this study. Correspondingly, among thirteen DEGs investigated, only five contribute to the increased butanol tolerance of the mutant B8 (Si et al., 2016). These results demonstrate that changes in the expression of a few genes may be valuable for functional improvements in butanol tolerance. Furthermore, the enhanced tolerance regulated by a TF may result from the combined actions of simultaneous changes in the expression levels of multiple genes. In addition, metabolic network and GO analyses of DEGs are necessary to identify potential mechanisms regulated by TFs, which would provide a basis for further elucidating the elaborate molecular interplay between a TF and its interacting proteins or responsive genes.

Conclusion

In summary, the rob gene-mutant strain DTrob was shown to exhibit higher butanol tolerance than BW25113, and 285 DEGs between the two strains were identified. Deletion of glgS and yibT and overexpression of fadB improved the butanol tolerance of E. coli. Furthermore, this is the first study to report that glgS can impact butanol tolerance directly. Analysis of the metabolic network of DEGs revealed that inactivation of Rob upregulated some genes that may be involved in the synergistic increase in acetyl-CoA production to promote cell growth under butanol stress. Thus, the results of this study provide a deeper understanding of the regulatory mechanism of Rob in the resistance of E. coli to butanol stress.

Data Availability Statement

The datasets generated for this study can be found in the Gene Expression Omnibus (GEO) of NCBI, GSE120032, https://www.ncbi.nlm.nih.gov/geo/query/acc.cgi?acc=GSE120032.

Author Contributions

YM conceived of the project, analyzed the data, and wrote the manuscript. YM, TX, and ZW designed the experiments and drafted the manuscript. TX, ZW, and DH performed the experiments. All authors analyzed the data, prepared the manuscript, and approved the final version.

Funding

This work was supported by the National Natural Science Foundation of China (NSFC-30900033), the Tianjin Science and Technology Council (18JCYBJC24200), and the Foundation (No. KF201815) of State Key Laboratory of Biobased Material and Green Papermaking, Qilu University of Technology, Shandong Academy of Sciences.

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.

Acknowledgments

We gratefully appreciate Dr. Yang Sheng (Chinese Acad Sci, Shanghai Inst Biol Sci, Inst Plant Physiol and Ecol, Key Lab Synthet Biol, Shanghai, China) for providing us with the plasmids and strains related to the CRISPR-Cas9 system.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fbioe.2020.524198/full#supplementary-material

References

Alexander, D. M., and St John, A. C. (1994). Characterization of the carbon starvation-inducible and stationary phase-inducible gene slp encoding an outer membrane lipoprotein in Escherichia coli. Mol. Microbiol. 11, 1059–1071. doi: 10.1111/j.1365-2958.1994.tb00383.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Aquino, P., Honda, B., Jaini, S., Lyubetskaya, A., Hosur, K., Chiu, J. G., et al. (2017). Coordinated regulation of acid resistance in Escherichia coli. BMC Syst. Biol. 11:1. doi: 10.1186/s12918-016-0376-y

PubMed Abstract | CrossRef Full Text | Google Scholar

Bennik, M. H. J., Pomposiello, P. J., Thorne, D. F., and Demple, B. (2000). Defining a rob regulon in Escherichia coli by using transposon mutagenesis. J. Bacteriol. 182, 3794–3801. doi: 10.1128/jb.182.13.3794-3801.2000

PubMed Abstract | CrossRef Full Text | Google Scholar

Black, P. N., Kianian, S. F., DiRusso, C. C., and Nunn, W. D. (1985). Long-chain fatty acid transport in Escherichia coli. Cloning Mapp. Exp. fadL Gene 260, 1780–1789.

Google Scholar

Brown, T. D. K., Jones-mortimer, M. C., and Kornberg, H. L. (1977). The enzymic interconversion of acetate and acetyl-coenzyme a in Escherichia coli. J. Gen. Microbiol. 102, 327–336. doi: 10.1099/00221287-102-2-327

PubMed Abstract | CrossRef Full Text | Google Scholar

Bui, L. M., Lee, J. Y., Geraldi, A., Rahman, Z., Lee, J. H., and Kim, S. C. (2015). Improved n-butanol tolerance in Escherichia coli by controlling membrane related functions. J. Biotechnol. 204, 33–44. doi: 10.1016/j.jbiotec.2015.03.025

PubMed Abstract | CrossRef Full Text | Google Scholar

Cai, L., Sutter, B. M., Li, B., and Tu, B. P. (2011). Acetyl-CoA induces cell growth and proliferation by promoting the acetylation of histones at growth genes. Mol. Cell 42, 426–437. doi: 10.1016/j.molcel.2011.05.004

PubMed Abstract | CrossRef Full Text | Google Scholar

Chin, W. C., Lin, K. H., Liu, C. C., Tsuge, K., and Huang, C. C. (2017). Improved n-butanol production via co-expression of membrane-targeted tilapia metallothionein and the clostridial metabolic pathway in Escherichia coli. BMC Biotechnol. 17:36. doi: 10.1186/s12896-017-0356-3

PubMed Abstract | CrossRef Full Text | Google Scholar

Dittrich, C. R., Bennett, G. N., and San, K.-Y. (2005). Characterization of the acetate-producing pathways in Escherichia coli. Biotechnol. Prog. 21, 1062–1067. doi: 10.1021/bp050073s

PubMed Abstract | CrossRef Full Text | Google Scholar

Dong, H. J., Zhao, C. H., Zhang, T. R., Zhu, H. W., Lin, Z., Tao, W. W., et al. (2017). A systematically chromosomally engineered Escherichia coli efficiently produces butanol. Metab. Eng. 44, 284–292. doi: 10.1016/j.ymben.2017.10.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Fujita, Y., Matsuoka, H., and Hirooka, K. (2007). Regulation of fatty acid metabolism in bacteria. Mol. Microbiol. 66, 829–839. doi: 10.1111/j.1365-2958.2007.05947.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Gao, C., Zhai, J., Dang, S., and Zheng, S. (2018). Analysis of alternative splicing in chicken embryo fibroblasts in response to reticuloendotheliosis virus infection. Avian Pathol. 47, 585–594. doi: 10.1080/03079457.2018.1511047

PubMed Abstract | CrossRef Full Text | Google Scholar

Gu, Y., Jiang, Y., Wu, H., Liu, X., Li, Z., Li, J., et al. (2011). Economical challenges to microbial producers of butanol: feedstock, butanol ratio and titer. Biotechnol. J. 6, 1348–1357. doi: 10.1002/biot.201100046

PubMed Abstract | CrossRef Full Text | Google Scholar

Guo, Y., Lu, B., Tang, H., Bi, D., Zhang, Z., Lin, L., et al. (2019). Tolerance against butanol stress by disrupting succinylglutamate desuccinylase in Escherichia coli. RSC Adv. 9, 11683–11695. doi: 10.1039/c8ra09711a

CrossRef Full Text | Google Scholar

Guzman, L. M., Belin, D., Carson, M. J., and Beckwith, J. (1995). Tight regulation, modulation, and high-level expression by vectors containing the arabinose PBAD promoter. J. Bacteriol. 177:4121. doi: 10.1128/jb.177.14.4121-4130.1995

PubMed Abstract | CrossRef Full Text | Google Scholar

He, X., Xue, T., Ma, Y., Zhang, J., Wang, Z., Hong, J., et al. (2019). Identification of functional butanol-tolerant genes from Escherichia coli mutants derived from error-prone PCR-based whole-genome shuffling. Biotechnol. Biofuels 12:73. doi: 10.1186/s13068-019-1405-z

PubMed Abstract | CrossRef Full Text | Google Scholar

Hommais, F., Krin, E., Coppee, J. Y., Lacroix, C., Yeramian, E., Danchin, A., et al. (2004). GadE (YhiE): a novel activator involved in the response to acid environment in Escherichia coli. Microbiology 150(Pt 1), 61–72. doi: 10.1099/mic.0.26659-0

PubMed Abstract | CrossRef Full Text | Google Scholar

Horinouchi, T., Maeda, T., and Furusawa, C. (2018). Understanding and engineering alcohol-tolerant bacteria using OMICS technology. World J. Microbiol. Biotechnol. 34:157. doi: 10.1007/s11274-018-2542-4

PubMed Abstract | CrossRef Full Text | Google Scholar

Isken, S., and de Bont, J. A. M. (1998). Bacteria tolerant to organic solvents. Extremophiles 2, 229–238. doi: 10.1007/s007920050065

PubMed Abstract | CrossRef Full Text | Google Scholar

Jang, Y.-S., and Lee, S. Y. (2015). Recent advances in biobutanol production. Ind. Biotechnol. 11, 316–321. doi: 10.1089/ind.2015.0023

CrossRef Full Text | Google Scholar

Jiang, Y., Chen, B., Duan, C. L., Sun, B. B., Yang, J. J., and Yang, S. (2015a). Multigene editing in the Escherichia coli genome via the CRISPR-Cas9 system. Appl. Environ. Microbiol. 81, 2506–2514. doi: 10.1128/aem.04023-14

PubMed Abstract | CrossRef Full Text | Google Scholar

Jiang, Y., Liu, J. L., Jiang, W. H., Yang, Y. L., and Yang, S. (2015b). Current status and prospects of industrial bio-production of n-butanol in China. Biotechnol. Adv. 33, 1493–1501. doi: 10.1016/j.biotechadv.2014.10.007

PubMed Abstract | CrossRef Full Text | Google Scholar

Jones, D. T., and Woods, D. R. (1986). Acetone-butanol fermentation revisited. Microbiol. Rev. 50, 484–524. doi: 10.1128/mmbr.50.4.484-524.1986

CrossRef Full Text | Google Scholar

Kaczmarzyk, D., Anfelt, J., Sarnegrim, A., and Hudson, E. P. (2014). Overexpression of sigma factor SigB improves temperature and butanol tolerance of Synechocystis sp PCC6803. J. Biotechnol. 182, 54–60. doi: 10.1016/j.jbiotec.2014.04.017

PubMed Abstract | CrossRef Full Text | Google Scholar

Kim, M., Hatt, J. K., Weigand, M. R., Krishnan, R., Pavlostathis, S. G., and Konstantinidis, K. T. (2018). Genomic and transcriptomic insights into how bacteria withstand high concentrations of benzalkonium chloride biocides. Appl. Environ. Microbiol. 84:e00197-18.

Google Scholar

Kim, Y. J., Im, S. Y., Lee, J. O., and Kim, O. B. (2016). Potential Swimming Motility Variation by AcrR in Escherichia coli. J. Microbiol. Biotechnol. 26, 1824–1828. doi: 10.4014/jmb.1607.07058

PubMed Abstract | CrossRef Full Text | Google Scholar

Krivoruchko, A., Zhang, Y., Siewers, V., Chen, Y., and Nielsen, J. (2015). Microbial acetyl-CoA metabolism and metabolic engineering. Metab. Eng. 28, 28–42. doi: 10.1016/j.ymben.2014.11.009

PubMed Abstract | CrossRef Full Text | Google Scholar

Li, H. G., Ma, X. X., Zhang, Q. H., Luo, W., Wu, Y. Q., and Li, X. H. (2016). Enhanced butanol production by solvent tolerance Clostridium acetobutylicum SE25 from cassava flour in a fibrous bed bioreactor. Bioresour. Technol. 221, 412–418. doi: 10.1016/j.biortech.2016.08.120

PubMed Abstract | CrossRef Full Text | Google Scholar

Livak, K. J., and Schmittgen, T. D. (2001). Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) Method. Methods 25, 402–408. doi: 10.1006/meth.2001.1262

PubMed Abstract | CrossRef Full Text | Google Scholar

Ma, Z., Gong, S., Richard, H., Tucker, D. L., Conway, T., and Foster, J. W. (2003). GadE (YhiE) activates glutamate decarboxylase-dependent acid resistance in Escherichia coli K-12. Mol. Microbiol. 49, 1309–1320. doi: 10.1046/j.1365-2958.2003.03633.x

PubMed Abstract | CrossRef Full Text | Google Scholar

Nakajima, H., Kobayashi, K., Kobayashi, M., Asako, H., and Aono, R. (1995). Overexpression of the robA gene increases organic solvent tolerance and multiple antibiotic and heavy metal ion resistance in Escherichia coli. Appl. Environ. Microbiol. 61:2302. doi: 10.1128/aem.61.6.2302-2307.1995

PubMed Abstract | CrossRef Full Text | Google Scholar

Nielsen, D. R., and Prather, K. J. (2009). In situ product recovery of n-butanol using polymeric resins. Biotechnol. Bioeng. 102, 811–821. doi: 10.1002/bit.22109

PubMed Abstract | CrossRef Full Text | Google Scholar

Peralta-Yahya, P. P., Zhang, F. Z., del Cardayre, S. B., and Keasling, J. D. (2012). Microbial engineering for the production of advanced biofuels. Nature 488, 320–328.

Google Scholar

Perez-Rueda, E., Tenorio-Salgado, S., Huerta-Saquero, A., Balderas-Martinez, Y. I., and Moreno-Hagelsieb, G. (2015). The functional landscape bound to the transcription factors of Escherichia coli K-12. Comput. Biol. Chem. 58, 93–103. doi: 10.1016/j.compbiolchem.2015.06.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Purvis, J. E., Yomano, L. P., and Ingram, L. O. (2005). Enhanced trehalose production improves growth of Escherichia coli under osmotic stress. Appl. Environ. Microbiol. 71, 3761–3769. doi: 10.1128/aem.71.7.3761-3769.2005

PubMed Abstract | CrossRef Full Text | Google Scholar

Qureshi, N., and Blaschek, H. P. (2000). Butanol production using Clostridium beijerinckii BA101 hyper-butanol producing mutant strain and recovery by pervaporation. Appl. Biochem. Biotechnol. 84-6, 225–235. doi: 10.1385/abab:84-86:1-9:225

CrossRef Full Text | Google Scholar

Rahimpour, M., Montero, M., Almagro, G., Viale, A. M., Sevilla, A., Canovas, M., et al. (2013). GlgS, described previously as a glycogen synthesis control protein, negatively regulates motility and biofilm formation in Escherichia coli. Biochem. J. 452, 559–573. doi: 10.1042/bj20130154

PubMed Abstract | CrossRef Full Text | Google Scholar

Ran, F. A., Hsu, P. D., Wright, J., Agarwala, V., Scott, D. A., and Zhang, F. (2013). Genome engineering using the CRISPR-Cas9 system. Nat. Protoc. 8, 2281–2308.

Google Scholar

Reyes, L. H., Almario, M. P., and Kao, K. C. (2011). Genomic library screens for genes involved in n-butanol tolerance in Escherichia coli. PLoS One 6:e17678. doi: 10.1371/journal.pone.0017678

PubMed Abstract | CrossRef Full Text | Google Scholar

Reyes, L. H., Almario, M. P., Winkler, J., Orozco, M. M., and Kao, K. C. (2012). Visualizing evolution in real time to determine the molecular mechanisms of n-butanol tolerance in Escherichia coli. Metab. Eng. 14, 579–590. doi: 10.1016/j.ymben.2012.05.002

PubMed Abstract | CrossRef Full Text | Google Scholar

Royce, L. A., Yoon, J. M., Chen, Y. X., Rickenbach, E., Shanks, J. V., and Jarboe, L. R. (2015). Evolution for exogenous octanoic acid tolerance improves carboxylic acid production and membrane integrity. Metab. Eng. 29, 180–188. doi: 10.1016/j.ymben.2015.03.014

PubMed Abstract | CrossRef Full Text | Google Scholar

Rutherford, B. J., Dahl, R. H., Price, R. E., Szmidt, H. L., Benke, P. I., Mukhopadhyay, A., et al. (2010). Functional genomic study of exogenous n-butanol stress in Escherichia coli. Appl. Environ. Microbiol. 76, 1935–1945. doi: 10.1128/AEM.02323-09

PubMed Abstract | CrossRef Full Text | Google Scholar

Sandoval, N. R., and Papoutsakis, E. T. (2016). Engineering membrane and cell-wall programs for tolerance to toxic chemicals: beyond solo genes. Curr. Opin. Microbiol. 33, 56–66. doi: 10.1016/j.mib.2016.06.005

PubMed Abstract | CrossRef Full Text | Google Scholar

Shen, C. R., Lan, E. I., Dekishima, Y., Baez, A., Cho, K. M., and Liao, J. C. (2011). Driving forces enable high-titer anaerobic 1-butanol synthesis in Escherichia coli. Appl. Environ. Microbiol. 77, 2905–2915. doi: 10.1128/aem.03034-10

PubMed Abstract | CrossRef Full Text | Google Scholar

Si, H. M., Zhang, F., Wu, A. N., Han, R. Z., Xu, G. C., and Ni, Y. (2016). DNA microarray of global transcription factor mutant reveals membrane-related proteins involved in n-butanol tolerance in Escherichia coli. Biotechnol. Biofuels 9:114.

Google Scholar

Sikkema, J., de Bont, J. A., and Poolman, B. (1995). Mechanisms of membrane toxicity of hydrocarbons. Microbiol. Rev. 59:201. doi: 10.1128/mmbr.59.2.201-222.1995

CrossRef Full Text | Google Scholar

Tucker, D. L., Tucker, N., Ma, Z., Foster, J. W., Miranda, R. L., Cohen, P. S., et al. (2003). Genes of the GadX-GadW regulon in Escherichia coli. J. Bacteriol. 185, 3190–3201. doi: 10.1128/jb.185.10.3190-3201.2003

PubMed Abstract | CrossRef Full Text | Google Scholar

van den Berg, B., Black, P. N., Clemons, W. M., and Rapoport, T. A. (2004). Crystal structure of the long-chain fatty acid transporter FadL. Science 304:1506. doi: 10.1126/science.1097524

PubMed Abstract | CrossRef Full Text | Google Scholar

Wang, X. J., Zhang, X., Yang, J. T., and Wang, Z. X. (2018). Effect on transcriptome and metabolome of stacked transgenic maize containing insecticidal cry and glyphosate tolerance epsps genes. Plant J. 93, 1007–1016. doi: 10.1111/tpj.13825

PubMed Abstract | CrossRef Full Text | Google Scholar

White, D. G., Goldman, J. D., Demple, B., and Levy, S. B. (1997). Role of the acrAB locus in organic solvent tolerance mediated by expression of marA, soxS, or robA in Escherichia coli. J. Bacteriol. 179:6122. doi: 10.1128/jb.179.19.6122-6126.1997

PubMed Abstract | CrossRef Full Text | Google Scholar

Zhang, S., Yang, W., Chen, H., Liu, B., Lin, B., and Tao, Y. (2019). Metabolic engineering for efficient supply of acetyl-CoA from different carbon sources in Escherichia coli. Microb. Cell Fact. 18:130. doi: 10.1186/s12934-019-1177-y

PubMed Abstract | CrossRef Full Text | Google Scholar