Genome‐wide transcriptomic and proteomic analyses of bollworm‐infested developing cotton bolls revealed the genes and pathways involved in the insect pest defence mechanism

Summary Cotton bollworm, Helicoverpa armigera, is a major insect pest that feeds on cotton bolls causing extensive damage leading to crop and productivity loss. In spite of such a major impact, cotton plant response to bollworm infection is yet to be witnessed. In this context, we have studied the genome‐wide response of cotton bolls infested with bollworm using transcriptomic and proteomic approaches. Further, we have validated this data using semi‐quantitative real‐time PCR. Comparative analyses have revealed that 39% of the transcriptome and 35% of the proteome were differentially regulated during bollworm infestation. Around 36% of significantly regulated transcripts and 45% of differentially expressed proteins were found to be involved in signalling followed by redox regulation. Further analysis showed that defence‐related stress hormones and their lipid precursors, transcription factors, signalling molecules, etc. were stimulated, whereas the growth‐related counterparts were suppressed during bollworm infestation. Around 26% of the significantly up‐regulated proteins were defence molecules, while >50% of the significantly down‐regulated were related to photosynthesis and growth. Interestingly, the biosynthesis genes for synergistically regulated jasmonate, ethylene and suppressors of the antagonistic factor salicylate were found to be up‐regulated, suggesting a choice among stress‐responsive phytohormone regulation. Manual curation of the enzymes and TFs highlighted the components of retrograde signalling pathways. Our data suggest that a selective regulatory mechanism directs the reallocation of metabolic resources favouring defence over growth under bollworm infestation and these insights could be exploited to develop bollworm‐resistant cotton varieties.


Introduction
Plants and insects have coexisted for about 350 million years leading to the evolution of both positive and negative interactions (Gatehouse, 2002). Positive interactions include insect-mediated pollination, seed dispersion, etc. that offer mutual benefit to the insect and the host, while negative interactions include insect predation (Gatehouse, 2002) that often causes detrimental effects to the host. In view of the long standing relationship, it is quite obvious that plants have evolved a diverse set of stressspecific constitutive and/or inducible defence mechanisms to resist and coexist with the insect pests (Gatehouse, 2002). The response pattern in both the mechanisms might either be activated locally at the infected site, systemically in the uninfected regions or by both the aforementioned through signalling molecules (Gatehouse, 2002). Signal perception and activation results in a vast cascade of events at cellular and molecular levels ultimately contributing to the defence mechanism. Specialized defence mechanisms that protect plants from insects include physical barriers such as cell wall and cuticle (Kempema et al., 2007), cellular processes including lignifications, cross-linking of cell wall components, release of volatile and nonvolatile metabolites (Kempema et al., 2007) and molecular processes like activation of defence-related genes and pathways (Ramirez et al., 2009;Ryan, 1990). In addition, hormones, transcription factors (TFs) and redox regulators play major role in stress response and defence signalling. Hormones are secondary signals that amplify primary elicitor signals during biotic stress (Yang et al., 1997). Salicylic acid (SA), ethylene (ET), jasmonic acid (JA) and systemin are the major phytohormones that are often quoted as stress-specific signalling molecules (Arimura et al., 2005;Loake and Grant, 2007;Sun et al., 2011;Yang et al., 1997). Moran and Thompson (2001) suggested that there is a complex crosstalk among hormonal pathways that control the plant responses to wounds, insect pest and pathogen attacks. In addition to hormones, pathogen elicitors also activate TFs that interacts with the pathogen-responsive cis elements present in the promoters of the defence-related genes. Even a single pathogen elicitor is capable of activating multiple TFs that can interact with the cis elements present within the same or different promoter regions ultimately leading to stimulation of vast set of defence-related genes and gene products (Yang et al., 1997). Oxidative burst is one of the major processes that occur during biotic stress condition (Lamb and Dixon, 1997). Maintenance of the redox balance within the plant cell plays a crucial role in modulating redox sensitive genes and proteins including many TFs (Torres, 2010). Release of reactive oxygen species (ROS) such as O À 2 and H 2 O 2 induces many defence-related events including cell wall reinforcement through lignifications and cross-linking of glycoproteins in the extracellular matrix, activation of defence-related genes, molecules, etc. (Jabs et al., 1997). The above-mentioned cascade of events demand and consume considerable amount of energy. During stress conditions, plant systems manage their biological energy and resources by either partitioning or favouring the molecular machineries towards defence and/or growth. As a result, most plants do survive insect predation; however, it is often accompanied by reduced growth and yield penalty depending on the site of insect predation.
Cotton (Gossypium spp.) is the leading contributor of natural fibre and is an important source of textile commodity, oil and protein meal (Han et al., 2004;Mei et al., 2004). Cotton bolls are crucial tissues that harbour lint (textile fibre) which is of huge economic value. The effect of insect pest infestation followed by secondary infection could lead up to 80% loss in cotton fibre production (Oerke, 2006). Around 1326 species of insects have been reported worldwide as cotton pests, and among them, bollworm (Helicoverpa armigera) is the major pest that directly feed and destroy the developing fibre tissue within the cotton bolls (Dua et al., 2006;Matthews, 1994). Biotic stress induced through insect pest attack regulates cellular events majorly driven by expression changes of genes and their associated pathways. Earlier efforts to understand plant diseases caused by insect, pathogen infestations have identified certain genes and pathways involved in the biotic stress tolerance in cotton (Artico et al., 2014;Dubey et al., 2013;Gao et al., 2013). Systems level analysis at the transcript and protein levels more often reveals the near-complete status of an organism subjected to stress or disease conditions (Komatsu et al., 2009;Srivastava et al., 2013). Such analyses are yet to be employed to understand molecular and cellular mechanisms operational during cotton plant and bollworm interactions. Further understanding of these interactions using high-throughput approaches might reveal stress-induced responses and endogenous resistance mechanisms operational in the host. In this context, we have made an attempt to understand the mechanisms adapted by cotton plant during bollworm attack using both transcriptomic and proteomic tools. As bolls are the target site for fibre synthesis as well as bollworm feeding, we have performed comparative analyses of developing cotton bolls subjected to bollworm infestation. Our comprehensive genome-wide analyses have revealed several new and interesting insights about cotton plant and bollworm interactions. Knowledge gained through this study could be further exploited to develop bollworm-resistant cotton varieties.

Differentially regulated transcripts and proteins in bollworm-infested cotton bolls
Cotton plants were grown in the field conditions following common agronomic practices. Only bolls that were infected by cotton bollworm (H. armigera) insect larvae were used for further analysis (Figure 1a,b). To study the effect of insect stress in developing boll tissue, microarray-based transcriptome profiling and two-dimensional gel electrophoresis (2D PAGE) followed by MALDI TOF/TOF-based proteome analyses were carried out at different developmental stages (Figure 1b,c). Labelled mRNA was hybridized to Affymetrix cotton GeneChip Genome array. Identified transcripts with a false discovery rate (FDR) adjusted P value ≤0.01 and fold change ≥3 were considered as differentially expressed transcripts (DETs) (Figure 1d; Table S1). In total, 8694 transcripts comprising 39% of the total transcripts present on the cotton GeneChip showed differential expression under bollworm infestation. Transcripts were annotated using the Arabidopsis TAIR protein database version 10 through BLASTX with E value cut-off ≤e-10. Identification of transcripts related to TFs, phytohormones and signal transduction were attained using Arabidopsis transcription factor and Arabidopsis hormone databases, respectively, as mentioned in Materials and methods. Classification of DETs under different functional categories was attained using MIPS functional catalogue. Gene expression patterns in response to bollworm infestation were classified using hierarchical clustering (  Tables S4  and S5).
Two-dimensional gel electrophoresis (2D PAGE)-based proteome analysis of bollworm-infested (biotic stress-induced BS) and noninfested (control, CN) cotton bolls showed an average of 393 reproducibly detected protein spots across developmental stages (0, 2 and 5 dpa) (Figures 3a-c and S1a-c). At least two independent replicate gels were generated per sample (Control, CN vs bollworm infested, BS). Only those spots that were reproducibly detected were further considered for comparative analysis. Quantification of the protein spots was attained using the per cent volume criterion that corresponds to the expression level of the detected spot regions. Protein spots showing fold change of AE1.5 with a P value <0.5 as mentioned in Materials and methods were considered as differentially expressed spots. Comparative analysis revealed that around 35% of the detected spots (137 spots) were differentially expressed (AE1.5-fold), and among them, 98 spots were identified using MALDI TOF/TOF (Table S14). Further Gene Ontology (GO)-based annotation and functional classification of the DEPs under various categories were attained using BLAST2GO platform version 2.7 (Figure 2e,f).
Comparative analysis of the transcriptome and proteome data sets highlighted 37 overlapping unique accessions that were quantified at both the transcript and protein levels (Table S15). Among them, only 10 genes (accessions) were found to have similar expression pattern in at least one of the developmental stages that were analysed in this study. Such poor corelation among the transcriptome and proteome data sets could either be attributed to the post-transcriptional regulation of the identified genes or be the experimental limitations associated with the transcript and protein turnover measurements.

Bollworm infestation induces early and consistent response in developing cotton bolls
Analysis of the transcriptome and proteome profile of the bolls that were subjected for only 8 h of bollworm infestation (0 dpa) showed 1352 (15.55%) DETs and 115 (36.37%) DEPs (Figure 1d, f). Majority of those transcripts (75.73%) and proteins (89.5%) were exclusively up-regulated in biotic stress-induced bolls (BS) at 0 dpa (Figure 1d,f and 4). Cluster analysis of the DETs and DEPs revealed that the gene expression pattern gradually changed upon boll development (Figure 5a,b). Briefly, 42.07% of 2-dpa, 72.85% of 5-dpa, 43.46% of 10-dpa transcripts were upregulated, while 54.05% of 2 dpa and only 25.75% of 5 dpa proteins were up-regulated in BS-induced bolls (Figure 1d,f). Analysis of the transcriptome data sets revealed an increase in the number of up-regulated transcripts at 5 and 10 dpa ( Figure 1d). However, such pattern was not observed at protein levels as the proteome profile showed a steady decline towards development (0-5 dpa, Figure 1f). In support of such drastic decline in transcript and protein populations, infested bolls showed compromised growth in terms of size accompanied by infectionrelated symptoms such as browning and rotting (Figure 1b). The discrepancies among transcriptome and proteome pattern on the one hand and growth of infested boll on the other suggest that the plant system's continuous effort to encounter pest attack at transcript level is somehow not been translated to the protein level. Nevertheless, the above-mentioned discrepancies can also be attributed to the limitations associated with transcript and protein quantitation procedures.

Defence signalling involves major reprogramming of metabolic and biosynthetic pathways
Gene Ontology-based functional annotation revealed that 62% of the DEPs were enzymes involved in regulating metabolic processes such as carbohydrate (10%), amino acid (19%) and lipid (18%) metabolisms ( Figure 2e). Cellular component-based classification showed that DEPs were distributed among plastid (20%), mitochondria (14%), nucleus (15%), extracellular (19%)  Table 1. Hierarchical clustering clearly showed that there were certain population of transcripts and proteins that were down-regulated upon development and certain others that were consistently upregulated across the developmental stages (Figure 5a,b). To explain the reprogramming pattern, we have focussed on the key metabolic pathways that showed major differences in terms of their expression exclusively during biotic stress conditions. Differentially expressed transcripts and DEPs related to trehalose, raffinose, malate, starch and cell wall metabolism majorly accounted for the carbohydrate metabolism (Table S9). Briefly, up-regulated transcripts related to trehalose phosphate synthase (TPS) and trehalose phosphatase (TPP) that are involved in trehalose biosynthesis (Tables 2 and S9) were identified in this study (Wingler, 2002). Two unique isoforms of galactinol synthases (GolS) involved in galactinol synthesis were found to be consistently up-regulated in our data set (Tables 2, S4, and S9). Further, our study also showed the up-regulated transcripts of malate synthase (MS) and down-regulated malate dehydrogenase (MDH) transcripts and proteins (Tables S9 and S14). Malate synthase is involved in the synthesis of malate, whereas MDH catalyses the oxidation of malate to pyruvate and CO 2 . The above-mentioned pattern in turn correlates with stress-responsive accumulation of trehalose, raffinose and malate in the infested bolls. In plants, synthesis and degradation of nonstructural carbohydrates (NSCs) such as starch plays a major role in the regulation of carbon source availability during growth and unfavourable conditions (Sulpice et al., 2009). Data curation revealed the down-regulated transcripts involved in starch biosynthesis such as starch synthase, ADP-glucose pyrophosphorylase (AGPase) and up-regulated enzymes of starch degradation pathway such as starch excess 1 and starch binding domaincontaining glycoside hydrolase (Yano et al., 2005) (Table 2, S1). Further analysis showed the down-regulated group of carbohydrate active enzymes (CAZymes) that catalyse cell wall metabolism-related processes such as loosening, elongation, grafting and maturation (Table S7). Analysis of nitrogen, amino acid and protein metabolism genes revealed that glutamine metabolism was up-regulated and ubiquitin cascade-related genes were differentially regulated. Pathway curation revealed that glutamine synthase isoform 1 (GS1), glutamate dehydrogenase (GDH) and nitrate reductase (NiR) involved in nitrogen metabolism were found to be up-regulated (Tables S4 and S14). In case of ubiquitin cascade, we observed that ubiquitin ligase (E3) was up-regulated, whereas ubiquitin conjugating enzyme (E2) was down-regulated in BS condition at transcript and protein levels (Tables S4, S11, and S14). This in turn suggests the onset of rate limiting pattern or controlled proteolysis through ubiquitin-mediated protein degradation during stress. In addition to the metabolic enzymes, we also observed up-regulated members of both sugar and amino acid transporters throughout the developmental stages suggesting active transportation processes (Table S12). Further, our study also revealed that fatty acid and lipid metabolism-related genes were differentially regulated (Table S10). Pathway mapping and curation suggested stimulation as well as repression of different lipid precursor pathways. Briefly, majority of the up-regulated genes were involved in glycolipid and phospholipid metabolism including enzymes related to a-linolenic acid metabolism that ultimately lead to JA biosynthesis (Tables S10 and S14). On the other hand, we observed that down-regulated genes such as 3oxo-5-alpha-steroid 4-dehydrogenase (DET 2) and 24-sterol Cmethyltransferase (SMT2-2) were involved in sterol biosynthesis. Sterols are membrane lipids that serve as precursor molecules for brassinosteroid (BR) biosynthesis (Table S10).

Repression of chloroplast, mitochondrial metabolism and cellular growth is evident during stress
In this study, the redundant expression pattern observed in majority of the metabolic pathways included both up-and downregulated genes. However, in case of photosynthesis, most of the DETs and proteins were found to be majorly down-regulated throughout the developmental stages. Briefly, genes related to photosystem I, II, ATP synthase, light harvesting complex, chlorophyll binding proteins, etc. were consistently downregulated under BS condition (Tables S8 and S14). Also genes encoding enzymes involved in carbon fixation including ribulose-1, 5-bisphosphosphate carboxylase/oxygenase (RuBisCO) and members of tricarboxylic acid cycle were found to be downregulated (Table S8). However, we observed transcripts related to pentatricopeptide repeat containing protein (GUN1) were upregulated (Tables S2 and S11; Figure 6). Further, transcripts related to mitochondria such as ATP synthase (mitochondrial), phosphate transporter, mitochondria-associated membrane glycoprotein (MAM33), etc. were found to be down-regulated (Tables S9 and S12), while transcripts related to mitochondrial alternative oxidase (AOX) involved in alternative respiratory pathway were found to be consistently up-regulated (Tables 2  and S1) (Vanlerberghe and McIntosh, 1997). In addition, cellular growth-related genes such as cytoskeleton proteins, annexins, profilins and expansins were found to be down-regulated upon development (Tables S7 and S14). Analysis of cell cycle and DNA replication-related genes revealed that members of cyclin-  Table S14. (c) 2D Spot profile of representative proteins showing differential expression under bollworm infestation in comparison with their respective control bolls during boll developmental stages (0, 2, 5 dpa). 2D PAGE proteome profile of infested and control bolls during developmental stages are presented in Figure S1.  (Table S7).
Calcium and redox signalling pathways are stimulated to regulate host response Calcium (Ca 2+ ) is often quoted as the second important messenger that activates signalling cascades in response to various stimuli including biotic stress (Sanders et al., 2002). Genes encoding Ca 2+ /calmodulin-binding proteins, calcineurin Blike proteins, Ca 2+ binding EF-hand family proteins, Ca 2+dependent ATPases were found to be up-regulated under BS conditions (Tables 1 and S14). Protein kinase cascades account for the major contributors of induced immunity in plants. In this study, we observed the up-regulation of two major stressactivated protein kinases such as calcium-dependent protein kinases and mitogen-activated protein kinases. Also, other kinases such as SNF1(sucrose nonfermenting)-related protein kinase, calcineurin B-like interacting protein kinase-6 (CIPK6), ankyrin protein kinases were found to be significantly upregulated in BS condition. The above-mentioned Ca 2+ signalling molecules can be broadly classified into nonenzymatic sensor proteins and enzymatic proteins. These proteins together constitute a network that not only maintains the intracellular calcium levels but are also involved in other signalling events including protein activation. Among the Ca 2+ -dependent enzymes, majority were observed to be kinases involved in phosphorylation reaction that determines the activity of substrate protein molecules. So we further analysed the data sets for kinase substrates and the pathways regulated by these kinases during stress. Interestingly, we observed CDPK-activated protein substrates such as NADPH oxidases and oxidoreductases that were up-regulated during later stages of development (10dpa) (Dubiella et al., 2013). Manual curation of data sets revealed that members of peroxidases, superoxide dismutases (SODs), glutathione S-transferases (GSTs), etc. were found to be differentially regulated at transcript and protein levels (Tables S13 and S14). These enzymes are actively involved in the ROS  Tables S2 and S3.   (Table S2). Further, around 8.1% (705) of the upregulated and 3.94% (343) of the down-regulated transcripts correspond to transcription factor families (TFs). Stress-responsive TFs belonging to WRKY, AP2-EREBP, NAC, bHLH, MYB, C2H2 and ethylene insensitive three families were found to be upregulated (Figure 6a). Further, a number of development-related TFs belonging to AUX/IAA, C2C2, GRAS and HB families were down-regulated during boll developmental stages (Figure 4a). Among the TFs, WRKY family accounted for about 13% of differentially expressed TFs in all of the analysed developmental stages (0-10 dpa). Following WRKY, AP2-EREBP, NAC and MYB TFs were found to be relatively more in the transcriptome data set. Among the stress-related TFs, AP2-EREBP plays a central role in the abscisic acid (ABA)-dependent stress signalling pathways. Literature survey and manual curation of the TFs composition and regulation pattern revealed that the TFs such as WRKY 40, NAC along with AOXs as mentioned elsewhere constitutively account for the components of mitochondrial retrograde signalling pathways (Sophia et al., 2013;Van Aken et al., 2013). Further analysis revealed the up-regulated members of WRKY33 and EIN3 that act as suppressors of SA.

Bollworm attack induces synthesis of synergistically regulated phytohormones
Phytohormone-related DETs constituted for about 64.5% of 0 dpa and <25% of 2, 5 and 10 dpa. Classification and annotation of phytohormone-related transcripts showed that around 24% of them were related to ABA followed by auxin (AUX/IAA), ethylene (ET), BR, SA, gibberellic acid (GA), JA and cytokinin. Manual curation revealed that transcripts corresponding to zeaxanthin epoxidase that catalyses the first step of ABA biosynthesis were found to be down-regulated (Table S3) (Marin et al., 1996). Transcripts corresponding to tryptophan aminotransferase (TAA1), aldehyde dehydrogenase of indole 3-pyruvic acid pathway (IPA) and cytochrome P450 enzyme-CYP79B2 of the indole-3-acetaldoxime pathway (IAOX) were found to be down-regulated. In addition, flavin monooxygenase (YUC), tryptophan decarboxylase of the indole-3-acetamide (IAM) pathway, transcripts related to nitrilases were found to be upregulated (Table S3). The above-mentioned pathways (IPA, IAOX and IAM) ultimately lead to auxin synthesis in plants. The aforementioned expression pattern suggests that auxin biosynthesis is partially inhibited during bollworm attack. Further, transcripts related to 1-aminocyclopropane-1-carboxylic acid (ACC) synthase and ACC oxidase involved in ethylene biosynthesis were found to be consistently up-regulated (Table S4). Transcripts related to sterol biosynthesis were found to be downregulated (Tables S3 and S10). Sterols serve as precursor for BR synthesis. In addition, transcripts related to 3-oxo-5-alpha-steroid 4-dehydrogenase and BR biosynthetic protein DWARF1 involved in BR biosynthesis also were found to be down-regulated. Isochorismatase hydrolase (ISH) was the only SA biosynthetic pathway-related enzyme that was found to be up-regulated in our data set (Table S1). However, pathway annotation revealed that ISH catalyses the conversion of isochorismic acid to 2, 3- dihydroxybenzoic acid (DHB) and this reaction does not necessarily lead to SA biosynthesis in plants ( Figure 6). Interestingly, most of the crucial enzymes involved in JA biosynthesis including phospholipase A, lipoxygenase, allene oxide synthase (AOS), allene oxide cyclase (AOC) and acyl-coA oxidase were found to be up-regulated throughout the developmental stages in BSinduced bolls at transcript and protein levels (Tables S10 and S14). Expression pattern of the hormonal biosynthesis-related genes suggests that biotic stress has induced JA and ET that are reported to act in a cooperative fashion ( Figure 6). Likewise, ABA being the major contributor for DETs showed downregulated pattern leading to its suppression during biotic stress. Down-regulation of BR and partial stimulation of auxin provides clues about the rate limiting pattern for allowing growth during stress.
Bona fide defence molecules and pathways are stimulated in response to biotic stress The ultimate output of metabolic reprogramming, transcription factor regulation, ROS, Ca 2+ signalling, hormonal biosynthesis, etc. leads to stimulation and synthesis of defence molecules and processes. In our data set, defence-related proteins such as the members of pathogenesis-related protein family (PR4, PR10, osmotin and thaumatin), chitinase, beta-glucanase and proteinase inhibitors were found to be up-regulated at transcript and protein levels (Figure 3c; Tables 2, S4, and S14). In addition, transcripts related to polyamine biosynthesis pathway such as SAM decarboxylase, spermidine synthase, arginine decarboxylase also were found to be up-regulated. Polyamines are cited as phytohormone like molecules that accumulate in response to stress and they also play major role in defence (Gill and Tuteja, 2010;Hussain et al., 2011;Waie and Rajam, 2003). In addition, ROS regulating enzymes such as peroxidases (cytosolic and extracellular), SODs, NADPH oxidase and oxidative stress-specific GSTs were found to be up-regulated (Table S13).

Comparative analysis highlights concordant and discordant members of transcript-protein pairs
Comparative analysis revealed 37 unique accessions that were commonly identified in both the transcriptome and proteome approaches (Table S15). Gene Ontology-based annotation and classification of these genes revealed that majority of them were involved in metabolic, cellular and biosynthetic processes including amino acid, nucleotide and osmolyte metabolism, stress and defence response and phytohormone biosynthesis ( Figure S3a-c). Manual curation of the data sets showed two distinct groups of transcript and protein pairs such as the genes with concordant (similar) expression patterns and the genes with discordant (dissimilar) expression patterns at transcript and protein levels (Table S15). Concordant members included genes such as chaperonin, actins, phosphoglycerate kinase, gibberellin oxidase, MDH and SODs that were found to be down-regulated, while defence and stress response-specific genes such as chitinase, PR protein, protease inhibitor and carbonic anhydrase were found to be up-regulated across developmental stages. Discordant members included ATP synthases, RuBisCO, glutamine synthase, proteasome subunits, etc. that showed poor corelation in their transcript and protein expression patterns.

Validation of microarray and proteome data by qRT-PCR analysis
To validate the data, semi-quantitative real-time PCR (qRT-PCR) analysis was performed on 42 selected differentially expressed genes (32 up-regulated and 10 down-regulated) during boll  developmental stages under cotton bollworm infestation (Figure S2). The results showed that the expression patterns of transcripts and proteins observed through microarray and proteome analyses were in parallel with those obtained by qRT-PCR ( Figure S2).

Discussion
In the current study, transcriptomic approach has resulted in the identification of relatively more number of differentially expressed genes as compared to the proteomic approach. Nevertheless, the expression pattern of crucial phytohormone biosynthesis genes including S-adenosylmethionine synthase and AOC, cytoskeleton proteins such as annexin isoforms and actins, cytosolic ascorbate peroxidase involved in redox regulation, signalling and stressspecific response genes such as calcium-binding protein and GS1 were exclusively identified at the protein levels. Transcriptomic approach revealed a vast set of TFs which are otherwise difficult to be identified at protein levels. The concordant pattern observed among the transcript-protein pairs such as chitinase, PR proteins, carbonic anhydrase and protease inhibitors adds significance and direct evidence for active defence signalling during bollworm infestation in developing cotton bolls. Also, discordance observed among transcript-protein pairs might reflect true biological discordance that could be attributed to post-transcriptional regulations, protein/transcript stability, miR-NAs, etc. and this needs to be investigated further. On the whole, our data suggest that employing two complementary approaches have increased the overall coverage of the differentially expressed genes that in turn has aided in filling crucial gaps in the abovementioned processes.

Host-pest interactions induce synthesis of additional metabolic resources ensuring survival
Transcriptomic and proteomic data obtained in this study revealed major changes in the carbohydrate metabolisms such as the significantly up-regulated genes encoding TPS, TPP and GolS involved in trehalose and raffinose biosynthesis (  and down-regulated genes involved in cell wall metabolisms (Table S7; Figure 6). Trehalose is a nonreducing disaccharide that plays a major role in the stabilization of proteins and molecular structures during stress (Garg et al., 2002). Likewise, raffinose belongs to a family of oligosaccharides that accumulates during stress and acts as osmoprotectants (Unda et al., 2012). Both trehalose and raffinose family of saccharides are referred to as compatible solutes responding to stress conditions in plants (Zhou et al., 2014). In addition to that, carbohydrate active and associated proteins like CAZymes, AGPs and FLAs were found to be differentially regulated in our data set (Table S7). AGPs are a heterogeneous class of abundant proteoglycans localized in both cell wall and cytosolic regions (Kumar et al., 2013). Reactive oxygen species molecules (H 2 O 2 ) released during stress conditions, cross-link the AGPs to the cell wall leading to rigidity and thereby render protection against pest and pathogen invasion. The rigidity caused in the cell wall matrix also acts as a negative regulator of cell growth (Cleland and Karlsnes, 1967;Gille et al., 2009;Sadava and Chrispeels, 1973). Analysis of the nitrogen and amino acid metabolism highlighted additional clues on BS regulation. Briefly, GS1 and glutamate-ammonia ligase were found to be consistently upregulated throughout the developmental stages (Table S14). In plants, glutamine (Gln) serves as the primary source for inorganic nitrogen, N (NO À 3 and NH þ 4 ) that gets subsequently utilized for biosynthesis of major amino acids like Glu, Asp and Asn. Among the genes involved in nitrogen/amino acid metabolism, glutamine synthetase (GS), glutamate synthase, GDH and NiR play primary roles in the assimilation of NH þ 4 . Two isoforms of GS are reported in plants among which GS1 is induced and the other isoform, GS2 is suppressed during pathogen attack (Pageau et al., 2006). Interestingly, our study revealed the consistent up-regulation of GS1, GDH and NiR under BS condition. Among the above mentioned, GS2 and NiR are involved in primary nitrogen assimilation, whereas GS1 and GDH are involved in organic nitrogen remobilization (Pageau et al., 2006). In addition, GS1 and GDH are also cited as senescence-related markers in plants (Pageau et al., 2006). Expression pattern observed in the current study in turn suggests that nitrogen remobilization and senescence leading to stress regulation is active, whereas signals related to primary nitrogen assimilation leading to growth are not evidenced.
Lipid metabolism serves as one of the major contributor for energy, membrane biogenesis, signalling molecules, etc. Our study revealed a bias in the stimulation of certain lipid metabolic pathways. Briefly, linolenic acid metabolic enzymes were up-regulated, whereas sterol biosynthesis genes were down-regulated. Linolenic acid and sterols are lipid molecules colocalized within the plastid and thylakoid membranes (Schwertner and Biale, 1973). Our data showed that the factors related to linolenic acid metabolism leading to biosynthesis of stress-responsive JA were induced, while the genes related to sterol biosynthesis that leads to growth-related Br synthesis were down-regulated. Interestingly, cellular component-based annotation of the above-mentioned pathways and processes revealed that growth and defence-related molecules are colocalized within the same compartment such as chloroplast. However, under biotic stress, only defence-related factors are positively regulated leaving behind the growth-related factors (Figure 6). Such a switch over in the regulation of lipid metabolism at subcellular level further ensures resistance during bollworm attack.

Diverse pathway regulation and association delineates bollworm infestation-specific signalling pattern
In addition to above-mentioned metabolic pathways, signalling molecules such as calcium, redox regulators, phytohormones, TFs and protein kinases that play independent role through diverse pathways were found to be differentially regulated in our data set (Figures S3 and S4). Briefly, up-regulated members of Calcium (Ca 2+ ) binding proteins and Ca 2+ transporting ATPases are involved in maintaining cytosolic calcium levels during stress conditions. Redox regulators and oxidative stress regulators such as SODs, GSTs, peroxidases, NADPH oxidases, respiratory burst oxidase homologs (Rbohs), DHARs, etc. were found to be temporally regulated in our data set (Figures 6 and S4; Table S13). Among them, NADPH oxidases, Rbohs and extracellular peroxidases are the major regulators of the primary apoplastic oxidative burst during insect attack (Torres, 2010). These enzymes catalyse reactions leading to ROS release which further stimulates downstream enzymes like SODs and GSTs localized at other cellular components. Superoxide dismutases and GSTs catalyse detoxification reactions, while PCBRs are involved in antioxidant synthesis ultimately protecting cells from oxidative stress (Niculaes et al., 2014). Temporal expression of the ROS scavengers indicates a compromised pattern executed by cotton bolls in response to bollworm attack. Phytohormones are secondary signals that often regulate development and stress conditions in plants. Our study showed the positive regulation of stress-related hormones such JA and ethylene accompanied by repression or down-regulation of growth-related Auxin, BR ( Figure 5). Interestingly, we did not find genes related to SA biosynthesis which is also a biotic stress-specific hormone; however, positively regulated suppressors for SA synthesis such as EIN3 and WRKY 33 have been evidenced in the current study (Table 2; Figure 6). Such observations in turn suggest that bollworm attack favours synergistic JA-ethylene synthesis over the antagonistic SA. Further, our study also revealed the upregulation of nuclear encoded plastid and mitochondrial components including TFs and enzymes. Among them, WRKY TF family is often linked with biotic stress response and pathogenassociated molecular pattern (PAMP) (Eulgem and Somssich, 2007;Rushton et al., 1996). In addition, WRKY TFs also regulate the expression of nuclear encoded mitochondrial proteins (Van Aken et al., 2013). For example; up-regulated members of WRKY40 identified in the current study are known to regulate the expression of AOX enzyme during stress conditions (Ivanova et al., 2014). These factors along with ROS molecules together constitute the plastid and mitochondrial signalling pathway components that regulates nuclear gene expression related to cell cycle and growth. Such pattern further reveals the retrograde trend operational during biotic stress conditions ( Figure 6).

Bollworm infestation induces major reallocation of metabolic resources favouring defence over growth
In response to insect pest attack, the host plant initiates several layers of defence including PAMP-triggered immunity (PTIs), effector-triggered immunity (Dangl and Jones, 2001). Following stress perception, a series of signal transduction events that include metabolic pathway regulation, defence molecule synthe-  (Figure S3). The above-mentioned coordinated events suggest an energy intensive mechanism that needs to be executed in order to exert the defence response. Knowledge gained through the current study highlights that the source for such additional energy could be attained by suppressing growth locally (boll tissue growth). In short, we observed a hierarchy of factors and processes that specifically suppresses growth-related events and stimulates defence-related processes. Expression trend of carbohydrate, amino acid and lipid metabolism also reveals the synthesis of stress response-related molecules such as trehalose, raffinose, linolenic acid and suppression of growth-related factors such as cell wall elongation enzymes, and sterols. Suppression of fundamental processes such as photosynthesis and cell cycle not only retards further growth but also regulates the amount of ROS released by them and also preserves considerable amount of energy that could be channelized for defence signalling. Likewise, defence response such as lignifications, cross-linking of proteoglycans to cell wall matrix offers cell wall rigidification on the one hand and negatively regulates cellular expansion on the other. All these factors together suggest a major reallocation of metabolic resources favouring defence over growth through selective regulation of specific pathways and processes during bollworm attack.

Conclusion
The present study delineates boll-specific endogenous defence mechanisms adapted by cotton plants under bollworm attack. The vast number of coordinated events including stimulations and repressions of major biological processes ultimately suggest major reallocation of metabolic resources that favours defence over growth in developing cotton bolls. Taking such insights into account, strategies targeting stimulation of multiple phytohormones and better sustainment of defence as well as growth signals could aid in developing resistant varieties against insect pest.

Plant material and biotic stress treatment
Cotton (Gossypium hirsutum cv. Bikaneri Narma) plants were grown at Agricultural Research Station, Dharwad farm, Dharwad, during 2012-2013 Kharif seasons following recommended agronomic practices. Two separate plots of the same genotype were maintained with a space of 90 cm between rows and 20 cm between plants. The plots were covered with nylon nets to protect from any external pest incidence. Plot designated as control (no infection from any class of insects including bollworms) and plot designated as infested (infested with H. armigera) were protected in early stage (45 days after sowing) from incidence of sucking pests by spraying recommended insecticides. During peak flowering stage (65-85 days after sowing), 2nd-3rd instar larvae of H. armigera, raised on bendi (Abelmoschus esculentus L. Syn. Hibiscus esculentus) fruits in the Entomology laboratory maintained at 25°C in 65%-70% relative humidity on a 14/10-h light/dark cycle, were released on buds of cotton on the day of pollination. The buds with larva were covered using paper bag with proper aeration to prevent larvae movement from the bud. Such a set-up ensured maximum damage of bolls by larva. The cotton bolls used as control(s) were also covered with paper bags with pores in order to prevent predation by insect pests and to ensure similar microenvironment as that of biotic stress-induced bolls. Samples were collected after 8 h of infection and labelled as 0 dpa, likewise samples collected after 2 and 5 days of insect infestation were labelled as 2 and 5 dpa, respectively. After 5 days of infestation, the insect was removed; the bolls were collected after 10 days of further growth and were labelled as 10 dpa. Harvested cotton boll samples were frozen immediately in liquid nitrogen and stored at À70°C until further use.

Total RNA isolation and Microarray experiments
Infected bolls (complete) from 0 and 2 dpa and only infested portion of 5 and 10 dpa boll samples along with their respective controls were used for RNA extraction. In order to minimize plant to plant and mode of infection variations, boll samples were collected and pooled from five independent plants and considered as one biological replicate. Total RNA isolation, analysis and quality check were performed as previously described. Affymetrix Cotton GeneChip Genome array (Affymetrix, Santa Clara, California) having 23 977 probe sets representing 21 854 cotton transcripts was used for transcriptome analysis (http://www.affymetrix.com/catalog/131430/AFFY/Cotton-Genome-Array#1_1). Three biological replicates were maintained to test the reproducibility and quality of the chip hybridization. Microarray hybridization, staining and washing procedures were carried out as described in the Affymetrix protocols with minor modifications (Padmalatha et al., 2012). The arrays were scanned with a GeneChip scanner 3000.

GeneChip data processing and analysis
After scanning of each array, DAT, CEL, CHP, XML and JPEG image files were generated using GeneChip Operating Software platform. The CEL files having estimated probe intensity values were analysed with GeneSpring GX-12.6 software (Agilent Technologies, Santa Clara, California) to get DETs. The robust multiarray average algorithm was used for the back ground correction; quantile normalization and median polished probe set summarization to generate single expression value for each probe set. Normalized expression values were log 2 -transformed, and differential expression analysis was performed using unpaired ttest. The P values were corrected by applying the FDR correction (Benjamini and Hochberg, 2000). Differentially expressed transcripts with FDR corrected P value ≤0.01 and fold change ≥3 were included for further data analysis. The hierarchical clustering was performed using complete linkage method with Euclidean distance based on log fold change data compared to control samples using Cluster 3.0 (Eisen et al., 1998) to display the expression pattern and tree diagram of DETs. The DETs were annotated using NetAffx annotation data for Cotton GeneChip (http://www.affymetrix.com, release 26).

Functional annotation of probe sets and pathways
To obtain functional annotation of transcripts, the consensus sequences of probe sets present in the Cotton GeneChip were mapped to the Arabidopsis TAIR protein database version 10 (http://www.arabidopsis.org) by BLASTX with E value cut-off ≤e-10. To identify the putative TFs and transcripts related to phytohormone biosynthesis and signal transduction pathways, the consensus sequences of all probe sets presented in cotton GeneChip were searched against the Arabidopsis transcription factor database (http://plntfdb.bio.uni-potsdam.de, version 3.0) and Arabidopsis hormone database (http://ahd.cbi.pku.edu.cn, version 2.0), respectively, by BLASTX with E value cut-off ≤e-10. Differentially expressed transcripts were grouped into functional categories based on MIPS functional catalogue (http:// mips.gsf.de/projects/funcat). Further, MapMan software version 3.5.0 (http://gabi.rzpd.de/projects/MapMan/) was used to visualize the expression of differentially regulated cotton transcripts onto metabolic pathways (Usadel et al., 2005). The microarray data are deposited in the Gene Expression Omnibus (GEO) database (http://www.ncbi.nlm.nih.gov/geo) at the NCBI under the series accession numbers GSE55511.

Two-dimensional gel electrophoresis (2D SDS-PAGE)based proteome analysis
The total protein from cotton bolls was isolated using phenol extraction method. The protein pellets were dissolved in 2D SDS-PAGE rehydration buffer (7 M urea, 2 M thiourea, 2% CHAPS, 0.5% ampholytes, 40 mM DTT), and an aliquot of protein sample from two independent replicates was subjected to two-dimensional gel electrophoresis (2D SDS-PAGE) as described previously (Kumar et al., 2013), briefly for the first-dimensional separation, the sample was loaded onto a 13-cm immobilized pH gradient (IPG) linear (pI 4-7) strips (GE Healthcare Life Sciences, U.S.A), and isoelectric focusing was performed according to manufacturer's instructions. Strips were then equilibrated, and second-dimensional separation was carried out on 12% SDS-polyacrylamide gel (13 cm, 1.5 mm). Gels were stained with Coomassie blue staining to visualize the protein spots and were stored in 1% acetic acid at 4°C until further use. Gels were scanned using GE Image scanner III (GE Healthcare Life Sciences, U.S.A) through Labscan software version 6.0.1 and analysed using Imagemaster 2D Platinum software version 6.0.1 (GE Healthcare Life Sciences, U.S.A). Protein spot detection parameters were set as: Smooth: 3, Minimum area: 11 and Saliency: 200. Detected protein spots were manually reevaluated to remove artefacts such as dust particles and streaks. Reproducibly detected protein spots were quantified using the per cent volume criterion. The relative volume corresponding to the detected spot region was considered to represent the expression level. Protein spots that showed normalized expression values of AE0.6-, 1.5-fold (biotic stress/control) were considered for statistical evaluation. To define the significant difference, P value <0.05 was set through Student's t-test and one-way ANOVA. Protein expression values within the above-mentioned thresholds were considered as differentially expressed.

Protein identification, annotation and classification from 2D gel spots
Differentially expressed protein spots were subjected to in-gel tryptic digestion followed by MALDIT TOF-based identification procedure as previously described (Kumar et al., 2013). Peptide MS/MS spectrum processing was achieved through Flexanalysis software version 3 and database search using Biotools software version 3.2. The database search parameters were set as described: fragment masses were searched in three independent databases: they were (i) NCBInr database (06/03/2010) containing 10 551 781 sequences (total) including 290 173 sequences from green plants (Viridiplantae), (ii) Gossypium raimondii protein database containing 40 976 sequences downloaded from Cot-tonGen website (ftp://ftp.bioinfo.wsu.edu/species/Gossyp-ium_raimondii/CGP-BGI_G.raimondii_Dgenome/genes/), (iii) Gossypium arboreum protein database containing 40 134 sequences downloaded from Cotton Genome Project website (ftp://cotton:cotton321$@public.genomics.org.cn/Ca_all_Ver-sion2.GENE.pep.gz) through mascot search engine, taxonomy was set as Viridiplantae, enzyme was set as trypsin, fixed modifications included carbamidomethylation of cysteine, variable modifications included oxidation of methionine, protein mass was unrestricted, missed cleavage was set to 1, MS tolerance of AE100 ppm and MS/MS tolerance of AE/À0.75 da. Only peptides with an individual ion score of >40 (P < 0.05) were considered for protein identification. Identified proteins were sequences that were exported into BLAST2GO platform version 2.7 (www.blast2go.com/b2ghome) to attain GO-based annotation, classification and pathway mapping (Conesa et al., 2005).

Transcriptome and proteome data set integration
Unique protein sequences corresponding to the differentially expressed proteins were subjected to GEO Nucleotide Translated BLAST: tblastn analysis to obtain accessions corresponding to Gossypium hirsutum (taxid: 3635) transcripts. The search parameters for tblastn were set as follows: database-GEO; organism: G. hirsutum (taxid: 3635). Transcript accessions with E value cutoff ≤e-10 and >70% sequence identity were considered as matched sequence.

The quantitative real-time PCR (qRT-PCR) analysis
The qRT-PCR analysis was performed on selected differentially expressed genes to validate the microarray and proteome expression data. RNA isolation followed by cDNA synthesis, qRT-PCR analysis and fold change calculations were performed as previously described (Padmalatha et al., 2012). The list of primers used in the current study is presented in Table S6. The GhPP2A1 gene (accession no: DT545658) from G. hirsutum was used as reference gene to normalize the expression values (Artico et al., 2010).

Supporting information
Additional Supporting information may be found in the online version of this article: Figure S1 2D-PAGE profile of control (uninfected) and bollworm infested cotton boll proteome during developmental stages. 2D-PAGE profiles of total proteins obtained from control and bollworm infested (Biotic stress) cotton bolls: (a) 0 dpa, (b) 2 dpa cotton bolls, (c) 5 dpa cotton bolls. Equal amount of total proteins (500 lg) were loaded onto 13 cm IPG strips, pI 4-7 and protein samples were resolved using 12% SDS-PAGE gels. Figure S2 Validation of microarray and proteome data using qRT-PCR during boll development stages (0, 2, 5 and 10 dpa) of cotton under biotic stress. Y-axis represents the log 2 fold change values at various stages in the biotic stress as compared to their respective stages in control. Figure S3 Gene ontology based classification of commonly identified genes in transcriptome and proteome datasets under biological process (a), cellular component (b) and molecular function (c) categories. Key events in the signal transduction pathway activated in response to biotic stress (d).

Figure S4
Overview of gene expression changes in developing cotton bolls infested with bollworm. Table S1 Differentially expressed transcripts during boll development stages (0, 2, 5 and 10 dpa) under biotic stress as compared to their respective control samples.  Table S15 List of commonly identified genes in the transcriptome and proteome data sets. Table S16 Cluster wise list of differentially expressed transcripts included in the Hierarchical Cluster Analysis depicted in Figure 5a. Table S17 Cluster wise list of differentially expressed proteins included in the Hierarchical Cluster Analysis depicted in Figure 4b.