Insight into buffalo (Bubalus bubalis) RIG1 and MDA5 receptors: a comparative study on dsRNA recognition and in-vitro antiviral response.

RIG1 and MDA5 have emerged as important intracellular innate pattern recognition receptors that recognize viral RNA and mediate cellular signals controlling Type I interferon (IFN-I) response. Buffalo RIG1 and MDA5 genes were investigated to understand the mechanism of receptor induced antiviral response. Sequence analysis revealed that RIG1 and MDA5 maintain a domain arrangement that is common in mammals. Critical binding site residues of the receptors are evolutionary conserved among mammals. Molecular dynamics simulations suggested that RIG1 and MDA5 follow a similar, if not identical, dsRNA binding pattern that has been previously reported in human. Moreover, binding free energy calculation revealed that MDA5 had a greater affinity towards dsRNA compared to RIG1. Constitutive expressions of RLR genes were ubiquitous in different tissues without being specific to immune organs. Poly I:C stimulation induced elevated expressions of IFN-β and IFN-stimulated genes (ISGs) through interferon regulatory factors (IRFs) mediated pathway in buffalo foetal fibroblast cells. The present study provides crucial insights into the structure and function of RIG1 and MDA5 receptors in buffalo.

Introduction

Establishment of an immune response begins with the discrimination between “self” and “non self” components within an organism. To recognize “non self” pathogen-associated molecular patterns (PAMPs), vertebrates have evolved a number of germ-line-encoded receptors called pattern recognition receptors (PRRs). Over the past decades, several types of PRRs have been identified, the recent being cytosolic retinoic acid inducible gene 1 (RIG1)-like helicases (RLHs) including RIG1 and melanoma differentiation-associated gene 5 (MDA5) [1]. RIG1 and MDA5 along with another gene LGP2 (Laboratory of Genetics and Physiology 2) constitutes cytosolic RIG1-like receptors (RLR) family, capable of sensing double stranded RNA (dsRNA) signatures of viruses. RIG1 and MDA5 share sequence homology and transduce downstream signal through a common adaptor called mitochondrial antiviral signaling protein (MAVS) that activates interferon (INF) signaling pathways [2], [3], [4]. RLRs differ from other classes of PRRs viz. Toll-like receptors (TLRs) and NOD-like receptors (NLRs), in that these receptors lack repetitive receptor elements such as leucine rich repeats (LRRs) as recognition platform of PAMPs. Instead, they recognize viral dsRNA with a central DExD/H box helicase domain which is characterized by seven conserved motifs (I/Walker A, Ia, II/Walker B, III IV, V and VI) [5], [6]. This conserved “helicase” core of RIG1 and MDA5 is flanked by two tandemly arranged N-terminal caspase activation and recruitment domains (CARDs), and a cysteine-rich C-terminal regulatory domain (CTD) [7]. LGP2 has similar domain architecture as RIG1 and MDA5, but it lacks N-terminal CARDs [8]. In normal cell, CTD plays a role in autoinhibition of RIG1 and MDA5 through intramolecular interactions that mask the CARD signaling [9], [10]. Study on crystal structures of CTD suggests that it plays a prominent role in high affinity binding and selectivity for dsRNA [8], [11]. However, MDA5’s CTD displays lower magnitude of affinity towards dsRNA than that of full-length MDA5 [12]. Thus, despite homology, they differ in molecular mechanisms for dsRNA recognition and play non-redundant function in antiviral immunity [13], [14]. MDA5 detects long duplex RNAs in the genome of double-stranded RNA (dsRNA) viruses or dsRNA replication intermediates of positive-strand viruses, e.g. Encephalomyocarditis and vesicular stomatitis viruses [15], [16]. In contrast, RIG1 detects the 5′ triphosphate group (5′ppp) and blunt end of short dsRNA or single strand RNA (ssRNA) hairpins, which are often present in a wide range of negative strand viruses, (e.g. influenza A virus) and in some positive and double strand RNA viruses [17], [18]. The role of LGP2 in antiviral response is less clear as compared to RIG1 and MDA5. Since LGP2 lacks N-terminal CARDs, it cannot mediate MAVS based signal transduction. Nonetheless, it is considered to act as a suppressor of RIG1 and MDA5 elicited signaling [19], [20]. Other studies, however, have shown LGP2 as a positive inducer of RIG1 and MDA5 amplifying their responses via a positive feedback mechanism [21], [22].

Buffalo (Bubalus bubalis) is an important milk and meat producing animal in many Mediterranean and East Asian countries including India. Besides, Indian buffaloes are well-known for their exceptional ability of disease resistance [23], [24]. However, the genetic basis of disease resistance is still elusive due to paucity of information on molecular mechanism(s) of host-pathogen interaction. Here, we studied the mechanism of viral RNA recognition by buffalo RIG1 and MDA5 receptors through molecular dynamics (MD) simulations. We investigated the expression profile of corresponding genes constitutively, and in response to a stimulus in order to understand the mechanism of RLR mediated antiviral response in buffalo.

Materials and Methods

Materials

Taq DNA polymerase, 10X buffer, dNTP were purchased from New England Biolab (MA, USA) and Fermentas Maxima SYBR Green qPCR Master Mix (2X) was obtained from Thermo Fisher Scientific Inc. (PA, USA). Unlabelled goat anti-bovine polyclonal antibodies against RIG1 and MDA5 were obtained from Santa Cruz Biotechnology Inc, (CA, USA) and horseradish peroxidase-labeled rabbit anti-goat antibody was obtained from Sigma-Aldrich (St. Louis, MO, USA). Poly (I:C)/LyoVec HMW (endotoxin level<0.001 EU/µg) was obtained from Invivogen (San Diego, CA). Dulbecco Modified Eagle’s medium (DMEM), Trypsin-EDTA, foetal bovine serum (FBS) and Epidermal Growth Factor (EGF) were obtained from Sigma-Aldrich (St. Louis, MO, USA). L-Glutamine (Glutamax 100x) was purchased from Invitrogen corp., (Carlsbad, CA, USA). Penicillin-G and streptomycin were obtained from Amresco (Solon, OH, USA). Filters (Millex GV. 0.22 µm) were purchased from Millipore Pvt. Ltd., (Billerica, MA, USA). All other reagents were of analytical grade.

PCR Amplification

All PCR amplifications were performed in 25 µl reaction volume. Each reaction contained 2.5 µl 10X buffer, 200 µM of dNTPs, 0.5 µl of each primers (10 pmol), 0.5 units of Taq DNA polymerase and nuclease free water to bring the total volume to 25 µl. Around 100 ng of cDNA was used as template. The details of primers for different fragments of RIG1 and MDA5 genes are summarized in Table 1. Thermal cycling was performed on an ABI Veriti PCR machine where a touchdown protocol with annealing 58°C was followed (Table S1 in File S1). The PCR products were resolved on a 1.5% agarose gel (Figure S1 in File S1). The PCR products were cloned on to pTZ57R/T vector, plasmids were screened using PCR, and positive plasmids were custom sequenced.

Table 1

Primers used for amplification of buffalo RIG1 and MDA5 genes.

Gene/Fragment	Primer	Sequences (5′–3′)	Length (bp)	Tm (°C)
RIG1Fragment1	Forward	CAGGCATGACGGCAGAGCAGCG	22	62
	Reverse	TTGATAATGAGGGCATCATTGTATTTCCGCA	31	59
RIG1Fragment2	Forward	GTGTTTCAGATGCCAGACAAAGAGGAAGA	29	60
	Reverse	TCTTCCTCTGCCTCTGGTCTGGATCAT	27	61
RIG1Fragment3	Forward	GTCGCCGATGAAGGCATTGACATTGC	26	61
	Reverse	CCTGAGCCCAAGGGGACATTTCTGC	25	63
MDA5 Fragment1	Forward	AGATGTCGTCGGACGGGAGTTCAACAGA	28	63
	Reverse	CTGAATCACTTCCCATGGTGCCTGAATC	28	61
MDA5 Fragment2	Forward	GTTTGGCAGAAGGAAGTGTCAGCTGCT	27	61
	Reverse	CCTCATTGTACTTCCTCAGATGTTCTGCAC	30	62
MDA5 Fragment3	Forward	AGCTGCAAAAGAAGGAAATCGCAAAGATCGTG	33	62
	Reverse	GAAATCTTGAATCAAATGTCAGTCCTCATCAC	32	59

Immunohistochemistry (IHC)

Spleen, tonsil and thyroid samples were collected from Municipal Corporation of Delhi (MCD) Slaughter House, Ghazipur, New Delhi, India with permission for research use. The samples were fixed with 10% neutral-buffered formalin and embedded in paraffin. For antigen retrieval, thin sections of tissues (∼3 µm) were deparaffinized with xylene and rehydrated by washing with ethanol (95%) and MiliQ water. Non-specific sites were blocked by overnight incubation of sections with a blocking reagent (1% BSA, 10% FBS-PBS) at 4°C. Tissue sections were then stained with either anti-RIG1 (1∶100 dilution) or anti-MDA5 (1∶150 dilution) for 1 h at room temperature. Slides were washed three times with Tris-buffered saline (pH 7.6) supplemented with 0.05% Tween-20. Thereafter, the slides were treated with horseradish peroxidase-labeled rabbit anti-goat antibody for 1 h and washed as described above. The slides were then mounted with DAB (3, 3′-diaminobenzidine) substrate chromogen solution for 5 minutes. Counter-staining was performed with Mayer’s haematoxylin. The slides were rinsed with water, dried and visualized.

Poly I:C Treatment of Foetal Fibroblasts

Gravid uteri (2–3 months of pregnancy) were obtained from abattoir house, New Delhi Municipal Corporation. For poly I:C treatment, buffalo foetal fibroblast cells at about 70% confluence were harvested and seeded at a concentration of 1×105cells/ml into four well culture plates (Text S1 in File S1). The cells were cultured in DMEM (10% FBS) for 24 h at 37°C. After 24 h, serum containing medium was decanted and cells were washed twice with Dulbecco’s phosphate buffer saline (DPBS) solution. Cells were then maintained with serum free DMEM medium at 37°C for 6 h. The medium was subsequently discarded and cells were maintained either with serum free DMEM medium (control group) or serum free DMEM medium containing 2 µg/ml of Poly (I:C)/LyoVec (treatment group). The experiment was performed with four replicates for each treatment or control group.

Real Time PCR (qRT-PCR)

Real time PCR primers were designed by aligning gene sequences of several mammals including cow, pig, mouse, buffalo, and human (Table 2). Equal amount of RNA (quantified by Qubit fluorometer, Invitrogen), wherever applicable, were used for cDNA preparation (Superscript III cDNA synthesis kit; Invitrogen). All qRT-PCR reactions were conducted on a Light Cycler 480 II Real-Time PCR machine (Roche Diagnostics, USA). Each reaction consisted of 2 µl cDNA template, 5 µl of 2X SYBR Green PCR Master Mix, 0.25 µl each of forward and reverse primers (10 pmol/µl) and nuclease free water for a final volume of 10 µl. Each sample was run in duplicate. Analysis of real-time PCR (qRT-PCR) was performed by delta-delta-Ct (ΔΔCt) method [25].

Table 2

Primers used for mRNA quantitation of different genes by Real time PCR (qRT-PCR).

Gene	Primer	Sequence (5′-3′)	Tm (°C)	Product Size (bp)
RIG1	Forward	CTTGCAAGAGGAATACCACTTAAACCCAGAGAC	63	150
	Reverse	TTCTGCCACGTCCAGTCAATATGCCAGGTTT	63
MDA5	Forward	TCTGCTTATCGCTACCACAGTGGCAGA	61	135
	Reverse	TGCTCTCATCAGCTCTGGCTCGACC	63
LGP2	Forward	CCACCACGTCAATGTGAACCCCAACTTC	63	176
	Reverse	TGAGCACTGGCAGCTTCACTGACTTGTAGAT	63
MAVS	Forward	CATCAGGAGCAAGACACAGAACTGGGCA	63	115
	Reverse	AACGGGCCAAGGGCTGGAAGGAGACAC	66
IRF3	Forward	GATACTGCCCTGGCTGATATCTCAGCTG	63	235
	Reverse	GTTCAGGGCAGACCGGAAATTCCTCTTC	63
IRF7	Forward	GCCTCCTGGAAAACCAACTTCCGCTG	63	340
	Reverse	GCAGATGGTCCTCCAAGCAGCTCTG	63
IFNβ	Forward	CACCACAGCTCTTTCCAGGAGCTACA	61	168
	Reverse	GGAACTGCTGTTCTTGCTTCATCTCCTC	61
ISG15	Forward	AGATCAATGTGCCTGCTTTCCAGCAGCG	63	162
	Reverse	GACCCTTGTCGTTCCTCACCAGGATG	61
ISG54	Forward	AAGTGCACGGCAATCATGAGTGAGACCA	61	94
	Reverse	CCTCTACCAAGTTCCAGGTGAAATGGC	61
ISG56	Forward	CTATGTGAAACACCTGAGAGGCCAGAATGA	62	132
	Reverse	CAGGCATAGTTGCCCCAGGTAACCAG	63

Molecular Dynamics Simulations

The initial structures of RIG1 and MDA5 were obtained using homology modeling (Text S2 in File S1). Molecular dynamics (MD) simulations of modeled RIG1 and MDA5 were carried out in GROMACS 4.5.5 [26] using Charmm27 all-atom force field. All the structures were hydrated using SPC216 water molecules in two separate cubic boxes with a distance of 10 Å between the proteins and box edges. Physiological ionic strengths (0.15 M) of counter ions were added to neutralize total charge of the simulation systems. The electroneutralized systems were subjected to steepest descent energy minimization to remove steric conflicts between atoms. The minimized models were then subjected to MD simulations using NPT ensemble with constraints on the backbone heavy atoms for 100 ps, keeping the number of particles (N), the system pressure (P), and the temperature (T) constant. The temperature of the systems was coupled to 300 K via Berendsen temperature coupling scheme, and the pressure of the systems were maintained at 1 bar using Parrinello–Rahman algorithm. Finally, the models of RIG1, and MDA5 were subjected to 2 ns MD simulations without any constraints on the backbone atoms. Trajectory and subsequent data analyses were performed using VMD 1.9.1 [27], Grace 5.1.23 (http://plasma-gate.weizmann.ac.il/Grace/), and Pymol1.3rc2 (The PyMOL Molecular Graphics System, Version 1.3 Schrödinger, LLC) software. Principal component analysis of MD trajectories was performed using ‘g_covar’ and ‘g_anaeig’ programs of GROMACS (Text S2 in File S1).

Binding Free Energy Calculation

The binding free energies (ΔGbind) of RIG1-dsRNA and MDA5-dsRNA were calculated with GMXAPBS tool [28]. The detailed procedure is outlined in Text S2 in File S1.

Statistical Analysis

All statistical analyses were carried out using SYSTAT 13.1 software (SYSTAT Software Inc.) Analysis of variance (ANOVA) was used to test between groups and hour intervals. Fischer’s restricted least significant differences criterion was used to maintain the a priori type I error rate of 0.05.

Results

Sequence Analysis of Buffalo RIG1 and MDA5 Genes

The amino acid sequences of buffalo RIG1 and MDA5 were translated from respective nucleotide sequences identified in this study (GenBank Accession ID: KF517376 and KF517377; Text S3 in File S1). The putative conserved domains and critical binding site residues within RIG1 and MDA5 were identified using DELTA BLAST and Conserved Domain Database (CDD) of NCBI. Comparison of amino acid conservation among different domains of buffalo RIG1 and MDA5 sequences revealed that CARD, helicase, and CTDs of these genes share 31.1, 34.9 and 21.5% homology, respectively (Figure 1). Putative conserved domains along with conservation of critical binding site residues among buffalo RIG1, MDA5 and different model organism has been shown in Figure 1A. From Figure 1A it is clear that critical binding site residues of both RIG1 and MDA5 are conserved across mammals. Moreover, the residues responsible for ATP, Mg++, and RNA binding are almost identical. Phylogenetic trees constructed by maximum likelihood method (MEGA5.2 software) [29] revealed a similar pattern of evolution for RIG1 and MDA5 in different model organisms except perrisodactyla (horse), where, RIG1 was clustered with cetartiodactyla and MDA5 with primates (Figure 1B).

Figure 1

Amino acid sequence analysis of RIG1 and MDA5 receptors.

(A) Sequence comparison of functional domains and binding residues of RIG1 and MDA5. Multiple sequence alignment was performed by MAFFT web server and binding site residues (highlighted) were identified by DELTA-BLAST and CDD database of NCBI. (B) The evolutionary history was inferred by using the Maximum Likelihood method of MEGA5 based on the JTT matrix-based model. The trees with the highest log likelihood are shown. Initial tree(s) for the heuristic search were obtained automatically by applying Neighbor-Join and BioNJ algorithms to a matrix of pairwise distances estimated using a JTT model, and then selecting the topology with superior log likelihood value. Sequences used for analysis are provided in Table S2 in File S1.

Structural Analysis of RIG1 and MDA5 Proteins

The modeled buffalo RIG1 and MDA5 exhibits characteristic tertiary folds containing HEL-1, HEL-2, HEL-2i, and the CTD (Figure S2 in File S1). The stereochemical qualities of the modeled proteins were satisfactory, indicating that the models could be used in simulation studies (Figure S3 in File S1; Table S3 in File S1).

MD simulations were performed on free and dsRNA bound RIG1 and MDA5 to analyze their structural organization and intermolecular interactions. Analysis of MD trajectories revealed that the backbone as well as the individual residues of RNA bound proteins fluctuate less compared to free proteins (Figure S4 in File S1). This observation prompted us to study the global motions of individual domains during MD simulation. Principal component analysis (PCA) revealed that CTD of free RIG1 shows maximum displacements (Figure 2A). However, CTD of RNA bound RIG1 was comparatively rigid. Further, it was found that HEL-2i domain of both free and RNA bound proteins display movements (Figure 2A). As can be seen from Figure 3B, all four domains of RIG1, i.e. HEL-1, HEL-2, HEL-2i, and CTD showed rotations and translations upon RNA binding. Similarly, free as well as RNA bound MDA5 showed noticeable motions of individual domain (Figure 2B). In addition, the V shaped linker of RNA bound proteins showed major displacements as compared to unbound proteins. Collectively, all four domains were found to be involved in dsRNA recognition.

Figure 2

Structural movements of RIG1 and MDA5 receptors during MD simulation.

Porcupine plots showing structural movements of (A) RIG1 and (B) MDA5 receptors during MD simulations. The structural movements are shown for both RNA-bound and free proteins. The docked dsRNA molecule has been colored in violet color. Arrow heads indicate direction of motion and length of the arrows specifies the extent of displacement. The figure also shows the schematic representation of intermolecular interactions (C) RIG1 and (D) MDA5 with dsRNA after MD simulations.

Figure 3

Constitutive mRNA expression levels of buffalo RLR (RIG1, MDA5 and LGP2) genes in different tissues.

Relative abundance of mRNA was measured by qRT-PCR using RPS-18 as housekeeping control.

The intermolecular contacts of buffalo RIG1 and MDA5 with dsRNA was identified and compared with those of human (Figure 2C–D; Table 3). We found that buffalo and human share a similar pattern of dsRNA recognition mechanism. In buffalo and human RIG1, positively charged residues were primarily involved in dsRNA recognition. Additionally, hydrophobic and polar residues were found to interact with dsRNA. In contrast, polar and hydrophobic residues of MDA5 were mostly involved in dsRNA binding before MD simulation. However, the number of Arginine around dsRNA increased considerably after MD simulation (Table 3). On the contrary, MDA5 displayed greater number of Lysine around dsRNA, suggesting positively charged residues are indispensible for viral RNA recognition by buffalo RIG1 and MDA5.

Table 3

Comparative intermolecular interaction of RIG1 and MDA5 (human and buffalo) with dsRNA.

RIG1						MDA5
Human		Buffalo				Human		Buffalo
RIG1	dsRNA	Before simulation		After simulation		MDA5	dsRNA	Before simulation		After simulation
		RIG1	dsRNA	RIG1	dsRNA			MDA5	dsRNA	MDA5	dsRNA
Ile300	C13	Val299	C13	Val299	C13	Lys365	C8	Val364	C8	Val364	C8
Gly326	G14	Gly325	G14	Gly325	G14	Val366	C8	Gly390	G9	Thr392	A11
Thr347	C13	Gln511	A8	Thr346	G14	Asn419	G9	Thr411	G9	Thr411	G9
Lys379	C6	Arg641	C10	Gln348	C13	Lys450	G8	Arg603	G3	Arg603	G3
Arg637	C10	Arg668	U12	Lys378	C6	Gln584	G2	Thr722	C5	Cys721	G4
Arg637	G11	Arg670	C10	His380	G5	Lys587	G2	Arg723	C5	Arg723	C5 C6
Gly663	G11	His834	C1	Glu514	A8	Arg605	G3	Arg723	C6	Arg723	C5 C6
Arg664	U12	Trp912	G2	Gln515	G9	Arg728	C5	Gly751	C6	Gly751	C6
His830	C1			Lys522	U7	Ala757	G7	Ala752	G7	Ala752	G7
Glu510	A8			Arg641	C10	Gln771	G7	Lys998	C6	Gln763	C8
				Arg668	U12	Asn812	C9	Val1000	G5	Arg838	U12
				Arg668	C13	His927	C3	His922	C3	Lys884	G7
				Arg670	C10	Lys1002	C6	Met921	C4	Lys920	G10
				Asn724	C6	Val1004	G5	Asn578	U12	Lys920	A11
				Cys754	A8			Gln913	C8	Lys981	C4
								Glu919	A11	Lys998	C6
										Gly449	G7

The overall binding free energy (ΔGbind) of MDA5−dsRNA complex was found to be greater than that of the RIG1−dsRNA complex (Table 4). The affinity of dsRNA towards the receptors differs only at the electrostatic (ΔGcoul) and polar contributions (ΔGpolar). The ΔGcoul and ΔGpolar of MDA5−dsRNA complex were calculated to be greater than those of RIG1−dsRNA complex. Hence, it appears that polar and electrostatic interactions play dominant role for dsRNA recognition by RIG1 and MDA5.

Table 4

MM/PBSA binding free energies (kJ/mol) of RIG1/dsRNA and MDA5/dsRNA complexes.

Energy component	RIG1/dsRNA (kJ/mol)	MDA5/dsRNA (kJ/mol)
ΔGcoul	−23858.9 (6.0)*	−26288.4 (5.1)
ΔGps	28789.7 (6.2)	32055.2 (5.6)
ΔGpolar	4930.8	5766.8
ΔGvdW	−3431.18 (1.2)	−3487.69 (1.8)
ΔGnps	−96.0808 (0.3)	−93.3481 (0.5)
ΔGnonpolar	−3527.26	−3581.04
ΔGbind	−1403.5 (3.7)	−2185.78 (3.5)

*Standard errors are indicated in parenthesis.

ΔGcoul = Electrostatic energy.

ΔGps = Polar solvation energy.

ΔGpolar = Polar contribution (ΔGcoul+ΔGps).

ΔGvdW = van der Waals energy.

ΔGnps = Nonpolar solvation energy.

ΔGnonpolar = Nonpolar contribution (ΔGnps+ΔGvdW).

ΔGbind = Overall binding energy.

Constitutive Expression of RLR Genes in Different Tissues

Basal expression of the three buffalo RLR genes (RIG1, MDA5, and LGP-2) was studied in liver, spleen, intestine, kidney, and PBMC using qRT-PCR with RPS-18 (40S ribosomal protein S18) as the internal reference gene. It was found that the RLR genes were ubiquitously expressed in every tissue tested (Figure 3). The expression levels of RLR genes were comparatively higher in kidney and intestine than in spleen and PBMC. The ubiquitous nature of RLR expression of was also reflected at protein level. Immunohistochemical localization of RLR proteins were detected in primary immune organs (spleen and tonsil) as well as in thyroid, kidney, which are not directly involved in primary immune response (Figure 4). Marked staining can be seen for RIG1 and MDA5 across the red pulp area of spleen, follicular epithelial cells and parafollicular cells of thyroid tissue, renal tubular epithelial cells, mucosal and sub-mucosal layers of intestine. A lower expression of RIG1 and MDA5 was observed in crypts of tonsil and in hepatocytes. Negative control sections with an isotypic IgG antibody were run in parallel; no staining was observed in any of these control sections. The primary antibodies used for RIG1 and MDA5 were specific, as revealed from western blot analysis (Figure S5 in File S1).

Figure 4

Constitutive protein expression levels of buffalo RIG1 and MDA5 in different tissues.

Imuunohistochemical localization of RIG1 and MDA5 were detected in thin tissue sections using specific antibodies with DAB (3,3′-Diaminobenzidine) as substrate. Counterstaining was performed with Mayer’s haematoxylin.

Responses of RLR and Associated Genes to poly I:C

Following poly I:C exposure, mRNA expression pattern of RLR genes, downstream adapter genes (MAVS, IRF-3 and IRF-7) and the effector genes (IFN-B and ISGs) were analyzed in foetal fibroblast cells. Fibroblast cells provide an obvious advantage for studying RLR mediated antiviral response, because TLR system for viral recognition is dispensable in this cell type [30], [31] Primary culture was chosen over established immortalized cell lines because the later may have multiple mutations that can affect the permissiveness and responsiveness to viral infections [32]. All experiments were carried out on fibroblasts with less than 10 passages.

Expression Patterns of RLR Genes

As shown in Figure 5, RLR genes showed no change in expression after 1 h post treatment, but RIG1 and MDA5 genes were up-regulated (p<0.01) after 4 h post treatment. Maximum expression of these genes was seen during 6 h post treatment. Afterwards, the expression of RIG1 was sustained at a higher level until 48 h. The level of MDA5 dropped transiently during 12 h, but increased again during subsequent intervals of experiment. The expression pattern of LGP2 followed a slightly different pattern (Figure 5C). The mRNA expression level of LGP2 was significantly higher (p<0.01) during 6 h and was maximum at 24 h after poly I:C treatment. It appeared that the amplitude of response was higher for MDA5, while RIG1 maintained a sustained response following poly I:C treatment.

Figure 5

Poly I:C induced expression levels of RLR genes in buffalo foetal fibroblast cells over different time intervals.

Relative fold change in mRNA of (A) RIG1 (B) MDA5 and (C) LGP2 over control has been shown.

Expression Patterns of MAVS, IRF-3, and IRF-7

The expression pattern of MAVS showed an unusual trend following poly I:C treatment (Figure 6A). The mRNA expression level was significantly lower (p<0.05) during 4 h post-treatment. The expression level increased during subsequent hours and was significantly higher (p<0.01) during 12–24 h post-treatment. Again during 48 h the expression dropped to basal level. The induction of both IRF-3 and IRF-7 following poly I:C treatment was evident in the present experiment (Figure 6B–C). The mRNA levels of both IRF-3 and IRF-7 were significantly higher during 6 h post treatment, though amplitude of response was much higher for IRF-7. The expression level showed a marked increase (p<0.01) during 24 h for both IRF-3 and IRF-7. In general, the expression pattern of IRF-3 and IRF-7 followed a similar trend until 24 h, after which the level of IRF-3 was constant until 48 h, but the level of IRF-7 dropped close to basal level.

Figure 6

Poly I:C induced expression levels of RLR adapter genes in buffalo foetal fibroblast cells over different time intervals.

Relative fold change in mRNA of (A) MAVS (B) IRF-3 and (C) IRF-7 over control has been shown.

Expression Patterns of IFN-β and ISGs

Poly I:C was found to induce IFN-β response in bovine foetal fibroblast cells. As shown in Figure 7A, IFN-β expression started increasing after 4 h post-treatment. Expression of IFN-β mRNA was highest during 6 h, higher level was maintained up to 12 h, and finally returned to basal level at 24 h post treatment. Induction of Type-I IFN response leads to transcriptional induction of >300 IFN-stimulated genes (ISGs) that have been shown to function as antiviral effectors. This was evident in the present experiment as well. As early as 4 h post-treatment, the transcript level of all ISGs (ISG-15, ISG-54 and ISG-56) increased (Figure 7B–D). At 6 h post treatment, a sharp increase in the expression level of all three ISGs was observed, the level of ISG-15 being the maximum. It was found that during subsequent time period the levels of all ISGs were sustained at higher concentration.

Figure 7

Poly I:C induced expression levels of IFN-β and ISGs in buffalo foetal fibroblast cells over different time intervals.

Relative fold change in mRNA of (A) IFN-β (B) ISG15 (C) ISG54 and (D) ISG56 over control has been shown.

Discussion

Over the last decade, studies have highlighted the biological role and importance of RLRs in antiviral immunity [7]. However, compared to other mammalian, avian, amphibian, and fish [33]–[35] species, bovine RLRs have been poorly studied. Here, we reported the evidence of functional RIG1 and MDA5 genes in buffalo genome. Domain analysis revealed that buffalo RIG1 and MDA5 receptors have common domain architecture as in other mammals. It is noteworthy that only MDA5, but no RIG1 ortholog, is present in all teleosts, the first chordates to be evolved. This indicates that MDA5 might precede RIG1 in the evolutionary history. It has been postulated that the functional domains present in MDA5 and RIG1 have evolved independently, and most probably, not by a simple gene-duplication event [36]. The helicase domains might have been evolved from a eukaryotic ancestor gene, common to DICER and DEAD/DEAH box RNA helicase-I of MDA5 and RIG1, but the emergence of distinct CARD domains might have happened by domain grafting [36], [37]. However, identification of RIG1-like genes in the genomes of sea anemone (Nematostella vectensis) and sea urchin [38] suggests that RIG1 might have emerged from invertebrates, but gradually got lost in some fish lineages [37]. While the origin of RLR genes remains inconclusive, our sequence analysis indicated that the DEAD/DEAH box RNA helicase (HEL-1 and HEL-2) of these two genes may share a common ancestor as revealed by the conservation of critical binding site residues viz. ATP binding site [G-(C/S)-G-K-T and D-Q-R-R], Mg++ binding site [D-E-C-H], and nucleotide binding site [K-T-R-(A/Q)-(P/H)-(I/V/S/N)-T-(S/T)-V]. Despite high sequence identity at RNA binding sites (E-Q/H-I/V-K) of RIG1 and MDA5, substantial differences were observed within HEL-2i of both the proteins. However, the residues of HEL-2i domain were conserved in both RIG1 and MDA5 across other mammalian species considered for this study. Our result suggests that irrespective of their point of phylogenetic diversion, critical binding site residues of RIG1 and MDA5 are extremely conserved. This indicates that RNA binding pattern, a critical step for inducing antiviral innate immune responses, must be similar throughout the mammalian species.

High resolution crystal structures of RIG1 and MDA5 with or without viral RNA have provided a better understanding of RNA-induced receptor activation [5], [10], [39]–[42]. However, these studies were carried out on human, mouse, or duck models, but information on RNA mediated receptor activation in bovine is sparse. Our simulation studies elucidated a network of interactions between dsRNA and four domains of buffalo RIG1 and MDA5. This mode of interaction was consistent with those observed in the crystallographic structures [5]–[6], [39]. In both the receptors, positively charged residues were involved in providing affinity to dsRNA, as observed by the binding free energy calculation. Analysis of binding free energy terms indicated that MDA5 binds dsRNA with a higher affinity than RIG1 does. We found the number of positively charged amino acids that interacted with dsRNA were higher in MDA5 than RIG1. This could be the reason of difference in net binding free energy and affinity of dsRNA-receptor complexes.

The constitutive expression levels of the RLR genes were not restricted to the immune organs. Instead, organs like kidney and intestine showed a higher level of constitutive expression. This ubiquitous expression pattern of RLRs has also been observed in human, mice, and fishes [43], [44]. Higher expression level in intestine is consistent with these receptors’ role in maintaining the intestinal epithelial barrier that prevents entry of many viruses. Taken together, these findings suggest that many cell-types have the potential for RLR mediated response against viral invasion.

The inductive expression patterns of RIG1 and MDA5 in response to poly I:C were different. The induction of MDA5 started early, amplitude of expression was higher, but expression level varied over time. However, RIG1 was induced later, but its expression level was sustained during different hours of experiment. The higher amplitude of MDA5 induction could be attributed to its greater affinity for dsRNA than RIG1, as revealed from our MD simulation studies. It seems that LGP2 negatively regulated both RIG1 and MDA5, as its expression pattern followed a reverse trend compared to RIG1 and MDA5. The role of LGP2 in RLR mediated antiviral responses remains controversial as it has been reported both as a suppressor [19], [20] as well as a positive inducer [21], [22] of RIG1 and MDA5 elicited signaling. LGP2 lacks N-terminal CARD regions, essential for MAVS mediated signaling, and keeps a check on RIG1 and MDA5 mediated overproduction of IFNs [45]. Our results suggested that while RIG1, MDA5, and LGP2 are simultaneously involved in viral RNA recognition and the role of LGP2 might be regulatory in nature.

After sensing dsRNA, RIG1 and MDA5 interact with CARD domain of MAVS, leading to downstream activation of IRF-3, IRF-7 and IFN-β production [2], [46]. However, we found the mRNA expression of MAVS was significantly (p<.05) downregulated during 4 h. RLR activation can induce a selective proteasomal degradation of MAVS, which is essential for downstream signaling [47]. Also, poly I:C can significantly downregulate MAVS expression in early hours of induction by rapidly decaying the MAVS mRNA [48]. Therefore, the initial downregulation of MAVS could be poly I:C mediated degradation of mRNA. Subsequent increase may be a consequence of cells’ demand of MAVS for RLR signaling. MAVS activation causes homo and hetero dimerization of IRF3 and IRF7, which induce type I IFNs. IRF3 is constitutively expressed and involved in immediate antiviral responses, while IRF7 is expressed in a limited fashion and is induced after IFNβ production is initiated [49]. Consistent with these findings, the present study also showed delayed induction of IRFs. Expression of IRF7 was temporal with high amplitude, but IRF-3 consistently showed a higher expression level after induction.

Following poly I:C treatment, buffalo foetal fibroblast cells triggered an elevated IFN-β expression and transcription of IFN-stimulated genes (ISGs). The expression pattern of IFN-β was consistent with previous studies [50]. It is worth noting that despite negligible expression level of IFN-β during 24–48 h, the expression levels of ISG15 and ISG54 continued to be high. This could be explained as activated homo- and heterodimers of IRF3 and IRF7 are able to bind IFN-stimulated response elements (ISRE) at the promoter regions of ISGs and can direct IFN-independent activation of ISGs [51]. Among ISGs, ISG56 can negatively regulate ISG15, ISG56, IRF3 and IFN-β by interacting with the adapter protein MITA to disrupt virus-induced IRF3 activation and cellular antiviral responses [52]. In agreement to this finding, the expression pattern of ISGs in our study suggests that ISG56 could be a negative regulator of ISG15 and ISG54, besides inhibiting the expression of IFN-β.

Conclusion

In conclusion, functional RIG1 and MDA5 genes of buffalo were found to follow similar pattern of domain organization as in other mammals. MD simulation studies suggested that despite striking homology at binding site residues, RIG1 and MDA5 receptors have different affinity towards dsRNA. Buffalo RLR genes are ubiquitously expressed in different tissues, irrespective of function of the tissue as immune organ. Poly I:C can induce an antiviral response in buffalo foetal fibroblast cells through RLR mediated pathway.

Supporting Text, Tables, and Figures. Text S1: Methodology for isolation and culture of buffalo foetal fibroblast cells. Text S2: Procedures for model construction and simulation analysis. Text S3: Deduced amino acid sequence of buffalo RIG1 and MDA5. Table S1: PCR cycling parameters used for the amplification of RLR genes. Table S2: Sequences (NCBI Accession numbers) used for phylogenetic analysis. Table S3: Comparative analysis of the stereo-chemical parameters of RIG1, MDA5 and their corresponding templates. Figure S1: Agarose gel electrophoresis of amplified products of different fragments of buffalo RIG1 and MDA5 genes. Figure S2: Overall tertiary structure of helicase and C-terminal domains of (A) RIG1 and (B) MDA5 proteins of buffalo. The helicases and C terminal domains of buffalo RIG1 and MDA5 were modeled based on the X-Ray crystallographic structures of human RIG1 and MDA5 (termed as templates). The pair-wise alignment of buffalo and human sequences showed 82 and 84 percent sequence identities for RIG1 and MDA5, respectively. Due to this striking homology between target and template proteins, the built models retained all the key structural features of the templates that included HEL-1, HEL-2, HEL-2i, and the C-terminal domains. Figure S3: Ramachandran plot of buffalo (A) RIG1 and (B) MDA5 receptors. Stereochemical qualities of modeled proteins were found to be highly comparable to those of templates, indicating quality of the models were reasonably good to carry out further studies. Figure S4: Stability parameters of RIG1 and MDA5 receptors as a function of simulation time. (A) RMSD of RIG1. (B) RMSD of MDA5. (C) RMSF of RIG1. (D) RMSF of MDA5. In each panel, red color indicates stability parameters of dsRNA bound receptors and black color denotes receptors without dsRNA. Figure S5: Western analysis of intestinal tissue to determine specificity and cross reactivity of primary antibodies to RIG1 and MDA5. Dilution of primary antibodies 1∶500; Dilution of secondary anti goat IgG 1∶50000. The chemiluminescence was detected on x-ray films using Immobilon Western Chemiluminescent HRP (Millipore Corporation, MA, USA) as substrate.

(PDF)

Click here for additional data file.