Genome-wide identification and characterization of long non-coding RNAs conferring resistance to Colletotrichum gloeosporioides in walnut (Juglans regia)

Walnut anthracnose caused by Colletotrichum gloeosporioides (Penz.) Penz. and Sacc. is an important walnut production problem in China. Although the long non-coding RNAs (lncRNAs) are important for plant disease resistance, the molecular mechanisms underlying resistance to C. gloeosporioides in walnut remain poorly understood. The anthracnose-resistant F26 fruits from the B26 clone and the anthracnose-susceptible F423 fruits from the 4–23 clone of walnut were used as the test materials. Specifically, we performed a comparative transcriptome analysis of F26 and F423 fruit bracts to identify differentially expressed LncRNAs (DELs) at five time-points (tissues at 0 hpi, pathological tissues at 24 hpi, 48 hpi, 72 hpi, and distal uninoculated tissues at 120 hpi). Compared with F423, a total of 14,525 DELs were identified, including 10,645 upregulated lncRNAs and 3846 downregulated lncRNAs in F26. The number of upregulated lncRNAs in F26 compared to in F423 was significantly higher at the early stages of C. gloeosporioides infection. A total of 5 modules related to disease resistance were screened by WGCNA and the target genes of lncRNAs were obtained. Bioinformatic analysis showed that the target genes of upregulated lncRNAs were enriched in immune-related processes during the infection of C. gloeosporioides, such as activation of innate immune response, defense response to bacterium, incompatible interaction and immune system process, and enriched in plant hormone signal transduction, phenylpropanoid biosynthesis and other pathways. And 124 known target genes for 96 hub lncRNAs were predicted, including 10 known resistance genes. The expression of 5 lncRNAs and 5 target genes was confirmed by qPCR, which was consistent with the RNA-seq data. The results of this study provide the basis for future functional characterizations of lncRNAs regarding the C. gloeosporioides resistance of walnut fruit bracts.


Background
Walnut (Juglans regia L.) is a diploid tree species (2n = 32), with approximately 667 Mb per 1C genome and an N50 size of 464,955 (based on a genome size of 606 Mbp) [1]. It is an ecologically important 'woody oil' tree species worldwide [2], and its kernel is a rich source of nutrients with health benefits for humans [3]. The peptides extracted from walnut seeds have antioxidant and anticancer activities and have the protective effects on the oxidative damage induced by H 2 O 2 [4]. Recent advances in biotechnology and genomics show potential to accelerate walnut breeding, such as gamma-irradiated pollen inducing haploid walnut plants [5], constructing the novel Axiom J. regia 700 K SNP array [6], and combining different assemblies to obtain the optimal version [7]. Walnut anthracnose caused by Colletotrichum gloeosporioides (Penz.) Penz. and Sacc can cause leaf scorch or defoliation and fruit gangrene, which is currently the disastrous disease in walnut production [8]. Due to the long incubation period of anthracnose, the concentrated onset time, and the strong outbreak, the use of chemical fungicides is still the main method of disease control [9]. The C. gloeosporioides lifestyle transitions associated with the infection of the host include the following three stages: attachment, biotrophy, and necrotrophy [10]. The pathogen of C. gloeosporioides in walnut overwinters in the diseased part with mycelium, and begins to move when the temperature reaches 11-15°C in the following spring [11]. Specifically, the formation of adherent cells is critical for fungal development during the C. gloeosporioides infection [12]. In a previous study, LAC2 was revealed to contribute to the formation of adherent cells to enhance the pathogenicity of C. gloeosporioides [13]. However, it is unclear how walnuts recognize and resist infections by C. gloeosporioides, and the regulatory network of hub and peripheral genes underlying the resistance of walnuts to C. gloeosporioides remains uncharacterized. Therefore, elucidating the molecular basis of this resistance mechanism is imperative for the breeding of walnut resistant to C. gloeosporioides [8,14,15].
Long non-coding RNA (lncRNA) is a type of RNA comprising 200-1,000,000 nt and structural characteristics similar to those of mRNA, but it does not encode a protein [16]. The lncRNAs were initially considered to be the transcription 'noise' of protein-coding genes, and were often ignored in transcriptome analyses [17]. However, the continuous development of sequencing technologies and transcriptome analyses has revealed that many lncRNAs in Arabidopsis thaliana [18], Triticum aestivum [19], Zea mays [20], and other plant species are related to stress responses, morphological development, and fruit maturation. For example, a heat-responsive lncRNA (TCONS_00048391) is an eTM for bra-miR164a and may be a competing endogenous RNA (ceRNA) for the target gene NAC1 (Bra030820), with effects on bra-miR164a expression in Chinese cabbage (Brassica rapa ssp. chinensis) [21]. Qin et al. confirmed that the DROUGHT INDUCED lncRNA regulates plant responses to abiotic stress by modulating the expression of a series of stress-responsive genes [22]. In A. thaliana, two lncRNAs, COOLAIR and COLDAIR, are associated with FLOWERING LOCUS C and play an crucial role in vernalization [23,24].
Many recent studies have proved that lncRNAs are important for plant-pathogen interactions. A role for nine hub lncRNAs and 12 target genes in the resistance of Paulownia tomentosa to witches'broom was uncovered via a high-throughput sequencing experiment, and their functions were analyzed with an RNA-lncRNA co-expression network model [25]. In tomato (Solanum lycopersicum), the lncRNA16397-GRX21 regulatory network reportedly decreases the reactive oxygen species content and cell membrane damage to enhance the resistance to P. infestans [26]. Moreover, the involvement of the WRKY1-lncRNA 33,732-RBOH module in regulating H 2 O 2 accumulation and resistance to P. infestans was determined based on a comparative transcriptome analysis [27]. In cotton (Gossypium spp.), a functional analysis demonstrated that a lack of two hub lncRNAs, GhlncNAT-ANX2 and GhlncNAT-RLP7, enhances seedling resistance to Verticillium dahliae and Botrytis cinerea, possibly because of the associated upregulated expression of LOX1 and LOX2 [28]. In wheat (Triticum aestivum L.), lncRNAs have a tissue-dependent expression pattern that can respond to powdery mildew infections and heat stress [29]. Additionally, four kinds of lncRNAs have important effects on Puccinia striiformis infections [30]. However, there are no reports regarding the role of lncRNAs in the walnut fruit resistance to anthracnose.
In this study, Illumina HiSeq 4000 sequencing was used to analyze the disease-resistant (F26) and susceptible (F423) fruit bracts at different C. gloeosporioides infection stages. The number and characteristics of lncRNAs were analyzed. Additionally, the hub lncRNAs related to disease resistance were screened and functionally analyzed to predict the role of lncRNAs in walnut fruit bract resistance to anthracnose. To the best of our knowledge, this is the first report on walnut lncRNAs and their biological functions related to fruit bract resistance to C. gloeosporioides. Our data may be a useful resource for clarifying the regulatory functions of lncRNAs influencing walnut fruit resistance to C. gloeosporioides.

Results
Symptoms and physiological changes of walnut fruit infected by C. gloeosporioide The resistant (F26) and susceptible (F423) fruit bracts were infected by C.gloeosporioide, the fruit bracts of F423 showed obvious symptoms at 48 hpi; the diseaseresistant fruit F26 at 72 hpi. The susceptible samples showed obvious C.gloeosporioide conidial at 120 hpi (Fig. 1a). During the infection, the activities of some enzymes and the content of hormones also changed correspondingly. Compared to the F423, the activities of chitinase, ROS-scavenging enzymes (catalase, CAT and superoxide dismutase, SOD) and the content of H 2 O 2 in F26 were higher (Fig. 1b-e). The content of salicylic acid (SA) and jasmonic acid (JA) in F26 was significantly higher than that in F423, and reached a peak at 72hpi after infection (Fig. 1f, g).

Whole genome identification of lncRNAs expressed in walnut fruit bracts
To identify lncRNAs expressed in walnut fruits in response to C. gloeosporioides, we constructed 20 cDNA libraries from the anthracnose-resistant and the anthracnosesusceptible walnut fruits at the following five infection stages: tissue at 0 hpi (hours post inoculation), infected tissue at 24, 48, and 72 hpi, and distal uninoculated tissue at 120 hpi (Additional file 1: Table S1). The libraries were sequenced with an Illumina HiSeq 4000 platform. A total of 265.4 Gb clean data were obtained, with an average of 13.27 Gb per library. Approximately 69.7% of the clean reads in all libraries were mapped to the walnut reference genome (Additional file 2: Table S2). The aligned transcripts were assembled, combined, and screened with the FEELnc software to obtain 22,336 lncRNAs (length ≥ 200 nt, ORF coverage < 50%, and potential coding score < 0.5), including 18,403 unknown lncRNAs (23.97%) and 3933 known lncRNAs (5.12%) (Fig. 2a,b). The principal component analyses (PCA) revealed that the results at same infection point were parallel (Fig. 2c).

Characterization of walnut fruit bract lncRNAs
A total of 58,369 mRNAs and 22,336 lncRNAs were obtained for the walnut fruit bracts (all samples combined) (Additional file 3: Table S3, Additional file 4: Table S4). The lncRNAs were characterized according to their locations relative to the partner RNA. A total of 40,429 (67.57%) lncRNAs were located in intergenic regions (i.e., only 32.43% genic lncRNAs). Additionally, 19,767 (48.89%) and 7302 (37.63%) of the intergenic lncRNAs and genic lncRNAs were located in the antisense strand, respectively ( Fig. 3a) (Additional file 5: Table S5). Most lncRNAs contained two or three exons, which differentiated them from mRNAs (Fig. 3c). Moreover, there was considerable diversity in the distribution of mRNA and lncRNA lengths (Fig. 3b). Furthermore, the expression level of most lncRNAs was significantly lower than that of mRNAs (Fig. 3d).

Differentially expressed lncRNAs at various infection stages
The lncRNAs that were differentially expressed between the disease-susceptible F423 fruits and the diseaseresistant F26 fruits at different C. gloeosporioides infection stages were analyzed.  Table S7). These results revealed the similarities in the expression of lncRNAs and mRNAs. And the number of upregulated lncRNAs and mRNAs in F26 compared to in F423 was significantly higher at the early stages of C. gloeosporioides infection. The highest r value (0.77) for F423 was calculated for the MEdarkseagreen module and F423_48hpi (Fig. 5b). And the MEorange module contains 76 lncRNAs and 227 mRNAs (Fig. 5c). These results suggested that lncRNAs are closely related to the disease resistance of walnut fruit bracts.

Enrichment analysis of genes co-expressed with lncRNAs
The GO and KEGG pathway databases were used to analyze the genes co-expressed with lncRNAs in each significant module and MEorange module. In the MEviolet module, a total of 208 GO terms were assigned, including 106, 8 and 94 GO terms in "biological process", "cellular component" and "molecular functions", respectively (Additional file 9:Table_S9). Among these enriched GO terms, most of them were newly predicted lncRNAs (red) and 58,369 protein-coding transcripts (blue). c Distribution of exon numbers in protein-coding genes (red) and lncRNA genes (blue). d Expression levels of protein-coding genes and lncRNA genes presented as log 10 (FPKM + 1) values related to biosynthesis and gene expression regulation, and the ones related to plant immunity were "response to stimulus"(GO:0050896) (187 genes) and "cellular response to stimulus"(GO:0051716) (114 genes) (Fig. 6a).

Network analysis of hub lncRNAs
The hub lncRNAs are important for regulating the whole network. Therefore, we screened the 96 hub lncRNAs and In the MEviolet module, the 25 known target genes for 15 hub lncRNAs were found to be involved in multiple functions (Fig. 7a), such as probable galacturonosyl transferase 10 and ultraviolet-B receptor UVR8-like. In addition, target genes encoding receptor-like serine/threonine-protein kinase NCRK (XM_018958556.1) and eukaryotic translation initiation factor 5A-2-like (XM_018994862.1) are known resistance genes (Fig. 8a). In the MElightyellow module, 16 hub lncRNAs were generated and their 22 known target genes were involved in many functions (Fig. 8b). And the target genes encoding G-type lectin S-receptor-like serine/ threonine-protein kinase LECRK1 (XM_018950446.1), probably inactive leucine-rich repeat receptor-like protein kinase At2g25790 (XM_018989953.1) and TMV resistance protein N-like (XM_018961957.1) were known resistance genes (Fig. 8b). In the MEbrown2 module, 24 hub lncRNAs and their 15 known target genes were generated (Fig. 8c), the target gene encoding probable LRR receptor-like serine/threonine-protein kinase At3g47570 (XM_ 018962714.1) was konwn resistance gene (Fig. 8c). In the MEwhite module, 23 hub lncRNAs were generated and their 38 known target genes were involved in many functions (Fig. 8d). The target genes encoding putative disease resistance protein At1g50180 (XM_018965430.1), probable LRR receptor-like serine/ threonine-protein kinase At1g63430 (XM_018973294.1) and L-type lectin-domain containing receptor kinase IV.2like (XM_018954279.1) were konwn resistance genes (Fig.  8d). In the MEorange module, 18 hub lncRNAs were generated and their 24 known target genes were involved in many functions (Fig. 8e). And the target gene encoding the inactive LRR receptor-like serine / threonine-protein kinase BIR2 (XM_018967526.1) was konwn resistance gene (Fig. 8e). All disease resistance genes in walnut are listed in Additional file 12: Table_S12. These results suggested that lncRNAs may participate in the resistance of Fig. 6 a Significantly over-represented GO terms in violet module for target genes. b Significantly over-represented GO terms in lightyellow module for target genes. c Significantly over-represented GO terms in brown2 module for target genes. d Significantly over-represented GO terms in white module for target genes. e Significantly over-represented GO terms in orange module for target genes walnut bracts to C. gloeosporioides by acting on their target genes. Based on the enrichment results of KEGG, we predicted the possible pathway of hub lncRNAs (Additional file 13: Table_S13). Most of the hub lncRNAs and its target genes in the five modules are enriched in the pathways of material metabolism and biosynthesis. In the white module, the function of hub lncRNA pathway map showed that cyclicnucleotide-gated channels and MPK4, the target genes of lncRNA MSTRG.94840.7,were upregulated at 72hpi, which were enriched in "plant pathogen interactions" pathway (Fig. 9a). The target genes (SAUR and ABF) of lncRNA103441.8 were involved in "plant hormone signal transduction" pathway,which may be related to plant immunity (Fig. 9b).

Discussion
In previous studies, lncRNAs were identified and analyzed in various biological processes important for seed development [31], photomorphogenesis [32], fruit development [33,34], and biotic and abiotic stress responses [22,35]. Additionally, there has been substantial research on the role of lncRNAs in plant-pathogen interactions. In A. thaliana, lncRNAs reportedly enhance the resistance to Pseudomonas syringae pv. tomato DC3000 by promoting PR1 expression [36]. In tomato, lncRNA23468 functions as a ceRNA that modulates NBS-LRR gene expression by mimicking the target of miR482b, thereby increasing the resistance to P. infestans [37]. Walnut anthracnose has been responsible for the premature fruit drop and yield losses that have adversely affected walnut production in China [13]. In this study, we investigated the role of lncRNAs in the resistance of walnut fruit bracts to anthracnose based on sequence analyses. Walnut anthracnose is caused by C. gloeosporioides, which completes its infection process as a hemibiotroph [10,38]. First, conidia germinate to generate appressoria, which produce invasion pegs that initiate the infection into susceptible plants. The primary mycelium produced in plant cells exists as a biotroph, after which the secondary mycelium produced in the infected site switches to necrotrophic growth [39,40]. We previously determined that the C. gloeosporioides life cycle in walnut tissue involves attachment at 24hpi, biotrophy at 48hpi, and necrotrophy at 72hpi (data unpublished). In this study, RNA-seq was performed to build the lncRNA and mRNA profiles of the walnut fruit bract tissue at 0 hpi, infected tissue at 24, 48, and 72 hpi, and distal uninoculated tissue at 120 hpi. A total of 58,369 mRNAs and 22,336 lncRNAs were identified, including 3933 known lncRNAs and 18,403 unknown lncRNAs. Consistent with the results of similar studies on other organisms, the identified putative lncRNA had fewer exons, shorter transcripts, and lower expression levels than protein-coding genes [41,42].
The release of walnut reference genome [1], enabled the study of walnut genetics at a genome-wide scale. Based on the reference genome, the whole-genome resequencing [43], the development of high-density genotyping tools [44], and the genetic dissection of important agronomical traits in walnut [45] have been completed. The development of bioinformatic analysis technology has enabled researchers to reveal that lncRNA functions and characteristics are far more complex than previously thought [16]. A recent comparative transcriptome analysis between wild-type and WRKY1-overexpressing tomato plants revealed 199 lncRNAs (DELs) and indicated that many of the lncRNA target genes that are likely affected by WRKY1 and associated with the resistance of tomato to P. infestans are involved in the response to biotic stimulus (GO:0009607) and plant-pathogen interaction (KO4626) [26]. In another recent study, 4594 putative lncRNAs were identified in comprehensive dynamic lncRNA expression networks under heat stress conditions. Co-expression networks revealing the interactions among the differentially expressed lncRNAs, mRNAs, and microRNAs indicated that several phytohormone pathways are associated with heat tolerance, including salicylic acid and brassinosteroid pathways [21]. In the current study, we obtained 10,645 upregulated lncRNAs and 15,247 upregulated mRNAs among the five comparisons (F26_0hpi vs F423_0hpi, F26_24 hpi vs F423_24 hpi, F26_48 hpi vs F423_48 hpi, F26_72 hpi vs F423_72 hpi, and F26_120 hpi vs F423_120 hpi). The number of up-regulated lncRNAs and mRNAs in the F26 vs F423 was significantly higher at the early stages of C. gloeosporioides infection. The functions of lncRNAs cannot currently be inferred from their sequence or structure, but lncRNAs can function in trans mode to target gene loci distant from where the lncRNAs are transcribed [46]. In F26, a total of 5 modules related to disease resistance were obtained by WGCNA during the infection of C. gloeosporioides. Many target genes of lncRNAs in these modules are enriched in plant immune related items and pathways, such as "activation of innate immune response", "activation of immune response" in MElightyellow module, "defense response to bacterium, incompatible interaction" in MEbrown2 module, "defense response" and "immune system process" in MEwhite module. These results suggest that these genes may play important roles in the process of resistance to C. gloeosporioides of walnut fruit bracts. Phytohormones are known to be important in the regulation of defense responses in plants [47][48][49]. Plants can exhibit systemic acquired resistance through the salicylic acid (SA) / jasmonic acid (JA)-mediated signaling network [50][51][52][53]. In our study, a total of 32 genes were identified in the significantly enriched KEGG pathway "Plant hormone signal transduction". Meanwhile, there are 3 and 5 genes enriched in "jasmonic acid mediated signaling pathway"and "response to jasmonic acid" respectively. We also showed that some genes were enriched in "auxin-activated signaling pathway" and "cellular response to auxin stimulus" at 24 hpi. Therefore, auxin may play a role in the resistance of walnut bracts to C. gloeosporioides. In addition, our result showed that the phenylpropanoid biosynthesis was one of the most significantly enriched pathways in the process of resistance to C. gloeosporioides of walnut fruit bracts. In this pathway, phenylalanine ammonium lyase (PAL) is the key regulatory enzyme in altering the biosynthesis and accumulation of flavonoids and lignin [54]. Lignin plays a structural role in the secondary cell walls formation [55], and flavonoids mediate plants against UV radiation and act as a visual signal for attracting pollinators [56,57]. In Caragana korshinskii, C. korshinskii adjusts its phenylpropanoid biosynthesis process to water-deficit environments and activates PAL by drought stress [58].
During long-term evolutionary interactions with plants, several pathogens successfully cause effectortriggered susceptibility response (ETS) by producing a number of effectors. Simultaneously, plants have evolved R genes that recognize these effectors and function through highly specific interactions between effectors and their corresponding nucleotide-binding site and leucine-rich repeat (NB-LRR) class receptors [59]. In tomato, lncRNA23468 reportedly increases the expression of the NBS-LRR target genes (encoding R proteins), resulting in enhanced resistance to P.infestans [37]. In the current study, we detected 10 R genes among the target genes of 96 hub lncRNAs. During the infection of C. gloeosporioides on the walnut fruit bracts, the results of Fig. 10 Validation of selected lncRNAs and mRNAs in a quantitative PCR assay. Blue and red represent the F423 and F26 samples, respectively. Expression data were normalized against the data for the18S rRNA housekeeping gene and are presented as themean ± standard error; *p < 0.05, **p < 0.01 ). These findings imply that lncRNAs may help mediate the disease resistance of walnut fruit bracts through the target R genes. The specific interaction between lncRNAs and R gene needs further verification. The expression levels of five hub lncRNAs (MSTRG13585, MSTRG11713, MSTRG152205, MSTRG112028, and MSTRG62751) and their target genes were further confirmed by qPCR, the results of which were consistent with the RNA-seq data. The data presented here provides researchers with the biological basis for future investigations of the mechanism underlying the disease resistance of walnut fruit bracts.

Conclusions
In this study we generated the expression profile of lncRNA in anthracnose-resistant F26 and anthracnosesusceptible F423 at five times. Compared with F423, a total of 14,525 DELs were identified, including 10,645 upregulated lncRNAs and 3846 downregulated lncRNAs in F26. Bioinformatic analysis showed that the target genes of upregulated lncRNAs were enriched in immune-related processes, plant hormone signal transduction, phenylpropanoid biosynthesis and other pathways during the infection of C. gloeosporioides. Hub lncRNAs with high connectivity to disease resistant genes were predicted. These results contribute to our understanding of the potential mechanism by which lncRNAs involved in C. gloeosporioides resistance and will facilitate the functional verification of the lncRNA in the future.

Plant materials and fungal isolates
The scions of walnut seedling tree B26 was provided by walnut specialized farmers' cooperative of Dongliugang  [8,14], and it was found that B26 clone was highly resistant to anthracnose in fruit bract, and the 4-23 clone was highly susceptible to anthracnose in fruit bract. The fruits of B26 clone (i.e., F26) and 4-23 clone (i.e., F423) were used as experimental materials. The voucher specimen of F26 and F423 had been deposited to our lab but not to any publicly available herbarium. We didn't use wild plants in this study and according to national and local legislation, no specific permission was required to collect these plants. C.gloeospor-ioidesm9 isolates (GenBank ID: GU597322) used in this study were maintained by our group.

Fungal pathogen inoculation of walnut fruits
Colletotrichum gloeosporioides was cultivated on potato dextrose agar medium for 5-7 days at 28°C.To prepare conidial suspensions, the colonies were washed with sterile distilled water containing 0.05% (v/v) Tween 80, passed through a filter (40-100 μm pores), quantified with a hemocytometer, and diluted with sterile distilled water to 105-106 conidia/ml [0.001% (v/v) Tween 80 final concentration]. Healthy fruits from the east-, south-, and west-facing parts of each tree were collected in mid-June and disinfected with 0.6% sodium hypochlorite and rinsed with sterile water. The punch inoculation of the detached walnut fruits was completed as previously described [8]. Based on anatomical changes to the infected walnut fruit bract, samples of the inoculation site were collected at 0, 24, 48, and 72 hpi, 0 hpi as a control. Additionally, distal uninoculated tissue was collected at 120 hpi. Take two independent samples as biological replicates at each infection time (Additional file 1: Table S1). All samples were flash-frozen in liquid nitrogen and stored at − 80°C until analyzed.

Read mapping and transcriptome assembly
The quality of the raw sequencing data was checked with FASTQC (http://www.bioinformatics.babraham.ac. uk/projects/fastqc/). Adapters and low-quality tags in the raw data were eliminated. Ribosomal RNA data were also removed. The remaining clean reads for the 20 cDNA libraries were combined and mapped to the J. regia genome sequence (https://www.ncbi.nlm. nih.gov/genome/? term = Juglans% 20 regia) with the HISAT program (version 0.11.5) (parameter setting: -rna-strandness RF) [60].

Identification of lncRNAs
To obtain the potential long non-coding RNAs, based on all the assembled transcripts, we have firstly excluded the known transcripts according to the class code "=". Then the remaining transcripts were used to remove the potential protein coding transcripts, miRNA-like, and other transcript types via blasting against the database of Rfam, Refseq, Uniprot, miRbase, and Pfam. Finally, the remaining transcripts were employed for coding potential prediction by using FEELnc tool. First, the FEELnc filter was used to remove short transcripts (default 200 nt) and assess single-exon transcripts [63]. The FEELnc codpot predictors were used to calculate a coding potential score. The assembled sequences were used for reconstructing the transcriptome. Finally, RNAs longer than 200 nt and derived from ≥2 exons, with an ORF coverage < 50% and a potential coding score < 0.5 were designated as lncRNAs [64].

Classification of lncRNAs
The lncRNAs were analyzed regarding their corresponding positions in the reference genome and the positional relationships between lncRNAs and partner RNAs based on 10,000-100,000 fragments. The lncRNAs were then divided into genic lncRNAs (overlapping partner RNAs) and intergenic lncRNAs (lincRNAs). The genic lncRNAs were further divided as overlapping, containing, or nested subtypes. Intergenic lncRNAs were divided as divergent, convergent, and same strand subtypes.

Analysis of differential expression patterns
Genes differentially expressed between the diseaseresistant and susceptible fruits at five infection stages were analyzed with DESeq2 (version 1.22.1) [65]. After assessing the significance of any differences, the genes with a p value ≤0.05 and a |log 2 foldchange| ≥ 1 were designated as differentially expressed genes. The principal component analyses (PCA) of F26 and F423 were constructed using the prcomp() function shipped with the base R installation. The PCA result was visualized using the ggplot2 package in R.

Quantitative real-time PCR
Total RNA samples extracted from walnut fruitsat individual infection stages were analyzed by qPCR. Briefly, first-strand cDNA was obtained with the TransScript One-Step gDNA Removal and cDNA Synthesis Super-Mixfor qPCR (Transgen, China). The lncRNA expression level was quantified with the TransStart Tip Green qPCR SuperMix (Transgen) and the CFX Connect Real-TimeSystem (Bio-Rad). The qPCR program was as follows: 95°C for 30 s; 40 cycles of 95°C for 10 s and 60°C for 30 s. For a melting curve analysis, the temperature was increased from 70°C to 95°C (0.5°C/5 s). All samples were analyzed in triplicate. The 18S rRNA gene was used as a housekeeping gene. The cycle threshold (Ct) 2 -ΔΔCt method (Software IQ5 2.0) [66] was used for the relative quantification of mRNAs. The primers used for RT-qPCR were designed with Beacon Designer 7 software and were synthesized by Sangon Biotech (Shanghai, China; Supplementary Table S14).

Prediction of lncRNA functions based on co-expression
Co-expression modules were generated with the WGCNA package (version 1.67) as previously described [67] (http:// lab.genetics.ucla.edu/horvath/Coexpression Network/). The lncRNAs and mRNAs that were not detected in at least one infection stage were not considered. In this analysis, the soft thresholding power was set to 12, after which the adjacency function was used to construct the adjacency matrix. A topological overlap measure map was constructed based on the adjacency matrix to calculate the similarity matrix of the lncRNA and mRNA expression between different nodes. The lncRNAs and mRNAs were hierarchically clustered based on the algorithm. To generate a number of clusters, modules were defined after eliminating or combining branches. The co-expression module dynamic shear tree parameters were determined as described by Gerttula [68]. The minimum number of genes was set to 30, the split sensitivity (deep Split) was set to 2, and the other settings were software default parameters. The module was related to the trait, and the correlation matrix between the module and the trait was calculated. The module with the highest correlation coefficient and the smallest p value was designated as the module most relevant to the trait. In this study, a significantly correlated module was identified based on a correlation coefficient (r) ≥ 0.8 [64] and p < 0.05. The co-expression networks of lncRNAs and hub lncRNAs in highly correlated modules were generated with the Cytoscape software (version 3.7.1) [69].

Functional enrichment analysis
The genes targeted by lncRNAs were functionally annotated based on the GO and KEGG pathway (http://www. genome.jp/kegg/) databases. The KOBAS program (version 2.0) was used to determine the significantly enriched KEGG pathways among the target genes [70]. According to the operation requirements of KOBAS 2.0, All data files were written with a parser. The gene-term mapping can be retrieved by parsing the raw data files for each pathway. The gene annotation and gene-ID relations were retrieved from KEGG Genes and BioMart. We mapped the genes in all databases to KEGG GENES and KEGG ORTHOLOGY (KO). The gene-pathway and is stored in our backend SQL relational database. The FASTA protein sequence files were preprocessed for BLAST [71].