Molecular Characterization and Phylogenetic Analysis of Bovine Ephemeral Fever Viruses in Khuzestan Province of Iran From 2018 to 2020

Bovine ephemeral fever (BEF) is an arthropod-borne viral disease caused by the BEF virus (BEFV). This single-stranded RNA virus that affects cattle and water buffalo is endemic in tropical and subtropical regions including Iran. While BEF is a major disease of cattle in Iran, information regarding its agent, molecular characterization, and circulating viruses are highly limited. The current study aimed to, rstly, determine the genetic and antigenic characteristics of BEFV strains in Khuzestan province in Southwest of Iran from 2018 to 2020 and, secondly, to compare them with strains obtained from other areas. had the highest identity and the evolutionary among the studied strains. Multiple sequence alignment and G3 antigenic sites showed that these neutralizing epitopes are highly conserved among the strains of the Middle East cluster; however, the strains previously identied in Iran differed in three amino acids placed in G1 and G2 epitopes.


Introduction
Bovine ephemeral fever (BEF), also called three-day sickness, is a viral disease that affects cattle and water buffalo. BEF virus (BEFV) is endemic in tropical and subtropical regions of Australia, Africa, and Asia, especially the Middle East countries [1,2]. The disease is marked by a bi-phasic or multiple-phasic fever, anorexia, muscle stiffness, salivation, depression, lameness, ocular and nasal discharge, cessation of rumination, and constipation [3]. BEF causes considerable economic losses, mainly because of its consequences and imposing international trade restrictions [4,5]. BEF, caused by an arthropod-borne virus, usually occurs in the summer and early autumn months [3]. Weather conditions plays a unique role in the spread of the BEF. Rainfall, prevailing wind patterns, streams, or other ground water sources have a signi cant effect on its occurrence, particularly in wet and dry seasons [6].
BEFV belongs to the genus Ephemerovirus, Rhabdoviridae family and has a negative-sense singlestranded RNA (ssRNA) [7]. The length of the viral genome is 14,900 nucleotides and encodes ve structural proteins of the matrix, nucleoprotein, phosphoprotein, surface glycoprotein, and RNA-dependent RNA polymerase. Surface glycoprotein (G) is a class I transmembrane glycoprotein that includes the major neutralizing epitopes and can induce protective immunity in cattle [8,9]. Four distinct antigenic sites (G1-G4) have been determined on the surface of this protein [9][10][11]. The viral capsid binds to the cell surface via G protein and penetrates into the cells through receptor-mediated endocytosis. Therefore, G protein plays an important role in the virus replication process as well as induction of immune response in the host body [6].
BEF disease is prevalent in Middle East countries and is considered an enzootic disease in some areas of this region [4,12]. The rst detailed report of BEF in the Middle East was described in Egypt in 1909, in which the disease spread from Aswan, Nile Valley to Cairo, and spread across the Delta to the coast [11].
Afterward, the disease was reported in the Palestine and Jordan Valley in 1931 [12,13], followed by Israel, Turkey, Iraq, Saudi Arabia, Kuwait, Yemen, Syria, and Iran [12,[14][15][16]. In Iran, the rst report of BEFV isolation dates back to the occurrence of the virus in the southern and eastern regions of the country in 1974 [14]. The disease is mostly sporadic and is prevalent in southern and western areas with tropical climates [14]. In Iran, the latest outbreak is related to 2020. The disease was appeared in Khuzestan in March as sporadic cases and then gradually appeared in the other provinces and soon reached to at least fourteen provinces (unpublished data). The peak of the disease was observed in July and August. While no precise estimate is available about losses caused by the outbreak; 70% of livestock in the infected areas showed clinical signs of disease (unpublished data). Three other outbreaks are also recorded in 2006 [17] and in autumn of 2012 and 2013 [14].
Although BFE is known to occur as a single serotype worldwide, some antigenic variations are detected at virus epitope sites [18]. Detection of these variations would be useful for identifying the epidemiology of BEFV, genetic evolution, and design broad-spectrum vaccines [18]. To determine the genetic and antigenic variations of BEFV strains in Iran from 2018 to 2020, we analyzed the nucleotide and deduced amino acid sequences of complete G gene of BEFV isolates and compared the sequences with those of the GenBank database.

Detection of BEFV
RT-nested PCR showed that all suspicious samples were positive for BEFV and a speci c product of 420 bp was found. To determine the genetic characteristics of BEFVs circulating in Iran during 2018 and 2020, the glycoprotein encoding gene of 4 BEFV strains (4 strains related to 2018 and 4 strains related to 2020) were sequenced. Because the G gene sequences of the prevalent strains in each year were the same, only the sequence of one strain in each year was submitted to the GenBank. The accession numbers of MZ51169 and MZ51168 were acquired for 2018 and 2020 Iranian BEFV strains, respectively. Identity and evolutionary distance analysis G gene nucleotide and deduced amino acid sequences of Iranian BEFV strains were compared with other BEFV sequences in the GenBank using multiple sequence alignment. The IR-2018 and IR-2020 strains showed the lowest identity with South African strains (86.8-87.1% nucleotide identity and 93-93.5% amino-acid identity). The IR-2018 strain shared the most nucleotide identity with the IR-2020 strain (98.4%), and Indian strains (98.2%), respectively, while concerning amino acid sequence, it showed the highest identity (99%) with the 2020 Turkish, 2006 Israeli and IR-2020 strains. The G gene sequence of the IR-2020 strain showed the highest nucleotide identity with strains identi ed in Turkey in 2020 (99.4-99.6% identity) as well as Indian strains (99.1% identity). Two BEFV strains of 2020 Iran and Turkey were quite similar in terms of amino acid sequence and were more than 99% similar to the Indian strains. Two investigated Iranian strains had 98.4% and 99% nucleotide and amino acid sequence identity, respectively (Table 1 and Supplementary Table 1).
The evolutionary distance analysis of Iranian strains showed similar results with identity. These strains had the highest and lowest evolutionary distances with the South African and the Middle East strains, respectively (Table 1).

Phylogenetic analysis
The phylogenetic tree was constructed with the G ectodomain encoding sequences of the viruses. As shown in Fig. 1, the BEFV strains formed four major clusters of East Asian, Australian, South Africa, and Middle East strains. The IR-2018 and IR-2020 BEFV strains were grouped in the Middle East cluster with the Turkish, Indian, and Israeli strains. The viruses related to 2020 outbreaks in Turkey and Iran were placed on the same branch, while BEFV strains of Indian-2018 and 2019, Turkish-2020, and Iranian-2018 and 2020 formed a sub-cluster.

Analysis of antigenic sites
No codon insertion, deletion, or frame shift was found in the glycoprotein sequences of IR-2018 and IR-2020 BEFV strains. All changes have resulted from single nucleotide substitutions. Multiple sequence alignments of G1, G2, and G3 antigenic sites showed that these neutralizing epitopes are highly conserved among the strains of the Middle East cluster, while it was found that the three substitutions at positions 218 (K to R), 223 (E to P), and 503 (K to T in IR-2013 strain and E to T in IR-2012 strain) in Iranian BEFV strains previously identi ed compared to present strains. Positions 218 and 503 are putative N-linked glycosylation sites. If the entire sequence of the P gene is considered, there are at least 22 amino acid substitutions between the strains previously identi ed in Iran and those present (Fig. 2). The Japanese YHL strain differed at 15 positions from IR-2018 and IR-2020 strains in which two substitutions of aa 170 (N to T) and aa 503 (K to T) were located in G2 and G1 antigenic sites, respectively.

Discussion
Although BEF is one of the most important diseases of cattle in Iran [14], little is available about its agent and characteristics of circulating viruses. Therefore, in the present study, we characterized the prevalent BEFV strains in Khuzestan province southwest of Iran in 2018 and 2020.
In 2018, the incidence of BEF was sporadic and no report of the disease is available in other regions of the country. Our surveillance program of the BEF showed the absence of the disease in 2019. However, many animals suspected of contracting the three-day sickness were reported in Khuzestan province and other parts of the country in 2020. The RT-nested PCR results con rmed the occurrence of the disease in Khuzestan province as well as other regions of the country. The BEF disease epizootic has been reported in Turkey simultaneous with an outbreak in Iran [19]. The phylogenetic analysis, identity matrix, and distance evolution con rmed the high genetic closeness of prevalent strains in Turkey and Iran in 2020.
According to the phylogenetic analysis of the G gene ectodomain region, the BEFV strains were categorized into four clusters of East Asia, Australia, Middle East, and South Africa. However, previous studies have grouped the BEF strains into the three clusters of East Asia, Middle East, and Australia [20][21][22][23][24]. According to our results and the study of Omar et al [25], the BEFV strains of South Africa are distinct from those of other regions of the world and fall into the new cluster of South Africa. These strains showed the most distance evolution and the lowest similarity with the identi ed strains in the present study. Israel, so they strongly suggested that these strains reached India from Israel and Turkey while maintaining their genomic sequence [28]. The ndings of the present study also con rm this view. These viruses probably traversed across the countries like Turkey and Israel reaching Iran, and then India. Iran and India do not have a common border, so it is unclear exactly how the virus has been transmitted from Iran to India. The similarity between BEFVs from Iran's and Turkey's 2020-outbreaks (99.4-99.6% nucleotide similarity) was higher than that of IR-2018 and IR-2020 BEFVs (98.4% nucleotide similarity).
The simultaneous outbreaks of BEFV in Iran and Turkey and high nucleotide and amino acid similarity indicate the common source of these viruses.
Comparison of G1, G2, and G3 antigenic sites showed that these neutralizing epitopes are highly conserved among the Middle East strains; however, the strains previously identi ed in Iran differed in three amino acids placed in G1 and G2 epitopes, which two of them were putative attachment site of oligosaccharides. Japanese YHL strain (vaccine strain used in Iran) differed from the Middle East BEFVs, especially IR-2018 and IR-2020 strains in G1 and G2 epitopes. The amino acid substituted in the G1 epitope is a putative N-linked glycosylation site; therefore, it may affect conformation and recognition of the neutralizing epitopes. Consistent with the differences in antigenic sites, Iranian and YHL BEFVs were placed as two distinct clusters. Due to these differences in antigenic structure, the use of heterologous strains as a vaccine may not provide full protective immunity [29]. According to the similarity of circulating strains in the Middle East, the use of these strains to develop vaccines can play an important role in preventing this disease.

Conclusion
Comparison of the nucleotide and deduced amino acid sequences of the BEFV G gene related to different countries showed that these viruses are geographically divided into four clusters: Middle East, East Asia, South Africa, and Australia. Based on the chronology and geographical area, the outbreaks of Turkey (2020), Iran (2018 and 2020), and India (2018 and 2019) are proposed to be related, which also suggests that these strains can be used to develop a BEFV vaccine.

Ethics approval
Written consent was obtained from all farm owners before entering the study. Sampling was performed by a specialist veterinarian to minimize pain and injury to the animal. In addition, the current study is approved by the Ethics Committee of Shahid Chamran University of Ahvaz. All methods were performed in accordance with the relevant guidelines.

Collecting samples
This study is conducted in Khuzestan province, southwest of Iran, with a tropical climate. Heparinized whole blood samples were collected from cattle with clinical symptoms existing in affected farms. The most important clinical symptoms in these animals were persistent fever, muscle stiffness, constipation, ocular and nasal discharge, and spontaneous recovery after three days. A total of 50 and 40 samples were taken from suspicious animals in 2018 and 2020, respectively. Blood samples were centrifuged at 3000 rpm for 15 minutes to separate buffy coats. The buffy coat layer was transferred to a new tube and then rinsed three times.

RNA extraction and RT-nested PCR
RNA was extracted from buffy coats using a viral RNA extraction kit (Bioneer; South Korea) according to the manufacturer's protocol. Then, extracted RNA samples were reverse-transcribed using a kit (Yekta Tajhiz; Iran) to obtain cDNA. Screening of the positive BEFV samples was performed based on the ampli cation of the partial G gene via a nested PCR [30,31]. Primer sequences are indicated in Table 2.
Sequencing the full-length G gene BEFV G gene was ampli ed as described by Hsieh et al with some modi cation [32]. A nested PCR was designed according to primers introduced by Hsieh et al. In the rst run, G1F and G4R primers were used to amplify an 1872 length fragment with the following protocol: 94 C for 2 min, 25 cycles of 94 C for 50 s, 50 C for 50 s, 72 C for 75 s, and a nal extension of 72 C for 5 min. In the second run, the PCR product was used as a template and G1-G4 fragments were ampli ed. PCR products were directly sequenced or subcloned into the pTG19-T vector (Vivantis; Malaysia) using standard techniques and then sequenced again. Sequencing was performed by Bioneer Company (South Korea) with the same PCR primers in two directions. Sequences were trimmed with BioEdit software version 7.0.4.1 (mbio, Inc, North Carolina, USA). Obtained sequences were submitted to the GenBank and are available under the accession numbers MZ51169 and MZ51168.

Phylogenetic analysis
Available global nucleotide sequences of the BEFV G gene were retrieved from the NCBI Genbank database. All retrieved sequences from different countries and those of Iran generated in this study were aligned with ClustalW [33]. A phylogenetic tree was constructed based on G ectodomain encoding sequences using Molecular Evolutionary Genetics Analysis the MegaX software package (10.0.4) [34]. Tree construction was performed using Maximum-likelihood (ML) method, General Time Reversible (GTR) Model as the nucleotide substitution model [35], and gamma-distributed 4 (G4) according to a recent report [27]. Sub-tree Pruning & Re-grafting (SPR) branch swapping were applied to the analysis. The reliability of the branching was evaluated by the bootstrap analysis with 1,000 replicates.

Distance estimation and sequence identity matrix
The pairwise distance was estimated using MegaX software. The maximum composite likelihood was used as a substitution model. A discrete Gamma distribution (G) was used to model evolutionary rate differences among sites. The number of bootstrap replication was considered 1000. The sequence identity matrix was calculated using BioEdit software.

List Of Abbreviations
Bovine ephemeral fever (BEF); BEF virus (BEFV); single-stranded RNA (ssRNA); Maximum-likelihood (ML); General Time Reversible (GTR); gamma-distributed 4 (G4); Sub-tree Pruning & Re-grafting (SPR) Declarations Ethics approval and consent to participate Written consent was obtained from all farm owners before entering the study. Sampling was performed by a specialist veterinarian to minimize pain and injury to the animal. In addition, the current study is approved by the Ethics Committee of Shahid Chamran University of Ahvaz. All methods were carried out in accordance with relevant guidelines and regulations.

Consent for publication
Not applicable.

Availability of data and materials
The datasets used and/or analyzed during the current study are available from the corresponding author on reasonable request. Sequence data of this project have been deposited in the GenBank of the National Center for Biotechnology Information (NCBI) under the accession number MZ51169 and MZ51168.   Figure 1 The phylogenetic tree of G gene ectodomain sequences constructed by the Maximum-likelihood using the MegaX program. BEFV strains were grouped into four clusters of I: East Asia, II: Middle East, III: Australia, and IV: South Africa. Sequencing viruses in this study are marked with a black square and previously studied viruses in Iran are marked with a black circle.