Genome-wide association study and transcriptomic analysis reveal the crucial role of sting1 in resistance to visceral white-nodules disease in Larimichthys polyactis

Abstract

Introduction: Larimichthys polyactis is a promising marine fishery species, but visceral white-nodules disease (VWND) caused by Pseudomonas plecoglossicida causes significant losses. However, genetic resistance mechanisms to VWND remain elusive in this species.

Methods: This study combined genome-wide association study (GWAS) and transcriptome analysis to unravel resistance loci and transcriptional regulation in L. polyactis.

Results: As a result, GWAS on 946 infected fish genotyped by 100 K lipid chips identified 22 suggestive significantly associated single-nucleotide polymorphisms (SNPs), annotated 60 candidate genes, where DNA-sensing pathway were enriched. RNA-seq on liver tissues of resistant, sensitive, and control groups found immune-related pathways enriched in the comparisons of RL vs CL and RL vs SL, and autophagy-related pathways enriched in the comparisons of SL vs CL and RL vs SL. Then, the integration of GWAS and transcriptome analysis identified seven key genes associated with resistance to VWND. Among the genes, the expression levels of mRNA for genes related to the cyclic GMP-AMP synthase-stimulator of interferon genes (STING) signaling pathway, as well as the protein levels of STING1, were significantly upregulated in RL. Collectively, integrating KEGG pathway analysis, gene and protein expression analysis revealed that the importance of STING1 for VWND resistance.

Discussion: These findings deepen the available knowledge on molecular mechanisms of host genetic resistance to VWND and provide an important foundation for the selection and breeding of VWND-resistant L. polyactis.

Keywords: GWAS; Larimichthys polyactis; STING1; transcriptome; visceral white-nodules disease.

Figure 1

(A) Frequency distribution and daily feeding rate change curves under Pseudomonas plecoglossicida infection in Larimichthys polyactis. (B) Survival curve of 946 L. polyactis after P. plecoglossicida infection. The evolution in survival rate showed a trend of firstly slowly decreasing, then rapidly declining, and finally slowly decreasing to zero. (C) The dissection symptoms of visceral white-nodules disease in L. polyactis. Red arrows indicate white bacterial nodules in the liver and spleen.

Figure 2

Manhattan plot and Quantile-Quantile (Q-Q) plot of genome-wide association study (GWAS) for resistance to VWND in L. polyactis. Manhattan plot of GWAS for (A) binary status (BS) and (C) survival time (ST). The X-axis represents the chromosomes, and the Y-axis shows the –log₁₀(p). The red solid line represents genome-wide significance association threshold. The dashed line represents the suggestive significance association. Q-Q plot of GWAS for (B) BS and (D) ST. The X-axis represents the expected –log₁₀(p), and the Y-axis represents the observed –log₁₀(p). A total of 13 SNPs of suggestive significance association level were distributed on Chromosome 4, 7, and 14 in the BS trait. A total of 15 suggestive significantly associated SNPs were located on Chromosomes 4 and 9 in the ST trait.

Figure 3

(A) Principal component analysis of the control group (CL), sensitive group (SL), and resistant group (RL). Orange, red, and blue points represent each fish of the CL, SL, and RL, respectively. Ellipses represnt 95% confidence ellipses. PC1: first principal components; PC2: second principal components. (B) Number distribution and expression trend of differentially expressed genes (DEGs) among CL, SL, and RL. Red bar: Significantly upregulated DEGs; Blue bar: Significantly downregulated DEGs. (C) Venn diagram showing the overlap of DEGs among SL vs CL, RL vs CL and RL vs SL.

Figure 4

Difference of transcriptomic expression patterns among CL, SL, and RL after infection. (A–C) Bubble plot of the top 20 significantly enriched KEGG pathways of 1657 DEGs overlapped between SL vs CL and RL vs SL, 857 DEGs overlapped between RL vs CL and RL vs SL and 466 DEGs overlapped among SL vs CL, RL vs CL, and RL vs SL, respectively.

Figure 5

(A) Venn diagram between the 9,638 DEGs of RNA-seq and the 60 potential candidate genes of GWAS. The overlap genes were identified as the final candidate genes. (B) Protein–protein interaction network of candidate genes. Node circles size and colour are positively correlated with abundance of interacting genes. Line segments thickness and colour are positively correlated with protein interactions. Triangle means transcription factor. (C) Heat map analysis of the seven key genes in livers of CL, SL, and RL, respectively. Genes in red are upregulated while those in blue are downregulated.

Figure 6

Quantitative real-time polymerase chain reaction validation of the GWAS and RNA-Seq results. (A) Expression of 20 potential candidate genes from GWAS at 96 h after infection in L. polyactis. ** means an extremely significant difference (P< 0.01) compared with the control group. (B) Relative expression values from quantitative real-time polymerase chain reaction results and transcripts per million values in the RNA-seq were compared. X-axis represents the seven key genes name, and Y-axis represents the expression levels (Log₂(fold change)) in RL relative to SL. Data were presented as means ± SEM.

Figure 7

Comparison of expression patterns of hub genes in the cGAS-STING signaling pathway among CL, SL, and RL following infection. (A) cgas, (B) sting1, (C) tbk1, (D) irf3, (E) irf7, (F) nfκb1, (G) p65, (H) nfκb2, (I) ifnc, (J) ifnd, (K) ifnh, (L) ifnγ, (M) il-6, (N) il-1β, and (O) tnfα. The expression levels of cgas, sting1, tbk1, irf3, irf7, nfκb2, ifnc, ifnd, ifnh, ifnγ, il-6, il-1β, and tnfα were significantly upregulated in the RL group; nfκb1 mRNA expression levels showed no significant difference among the three groups; p65 mRNA expression levels was significantly downregulated in the RL group. Data were presented as means ± SEM. ** means an extremely significant difference (P< 0.01) compared with the control group or the SL group.

Figure 8

Comparison of protein expression level of STING1 among CL, SL, and RL following infection. Data were presented as means ± SEM. ** means an extremely significant difference (P< 0.01) compared with the control group or the SL group. The protein expression levels of STING1 in the RL group were significantly higher than those in the CL and SL groups.

Figure 9

The regulatory network of the canonical cGAS-STING signaling pathway.