All of the animal experiments were performed in the Jiangsu Academy of Agricultural Sciences (JAAS) and were reviewed and approved by the Committee on the Ethics of Animal Experiments of JAAS. All of the experimental procedures conformed to the guidelines of Jiangsu Province Animal Regulations (Government Decree No. 45) in accordance with international law.

M. synoviae, Mycoplasma gallisepticum, Mycoplasma gallinarum, Mycoplasma iowae, and Mycoplasma meleagridis strains were cultured in modified Frey's medium at 37°C with 5% CO₂ as described previously (24). Staphylococcus aureus strains were cultured in Todd-Hewitt medium at 37°C. Avian pathogenic Escherichia coli (APEC), Salmonella enterica subsp. enterica serovar Gallinarum biovar Pullorum (S. Pullorum), Salmonella enterica subsp. enterica serovar Gallinarum biovar Gallinarum (S. Gallinarum), and Salmonella enterica subsp. enterica serovar Typhimurium (S. Typhimurium) strains were cultured at 37°C in Luria-Bertani medium. The strains used in this study were passed fewer than five times in vitro.

A total of 46 strains from 25 major animal Mycoplasma species/subspecies that had complete genome sequences and annotation information published in GenBank were selected and used for comparative genomic analysis. Among these 46 strains, there were 18 strains of nine avian Mycoplasma species, including six strains of M. synoviae. General characteristics of the complete genomes of the 46 strains of the 25 species/subspecies of Mycoplasma are listed in Supplementary Table 1.

The translated amino acid sequence sets corresponding to all protein-coding genes were obtained from the RefSeq annotation entries of the 46 strains in GenBank and used for core/pan gene analysis. The core/pan genes of the 46 Mycoplasma strains were clustered by CD-HIT (v 4.8.1) rapid clustering of similar proteins with a threshold of 50% pairwise identity and 70% length difference cut-off in amino acid sequence (25). The evolutionary relationships among the 46 strains of the 25 species/subspecies of Mycoplasma were characterized via analysis of 16S rRNA gene sequences and all single-copy core gene nucleotide sequences. All of the sequences mentioned above were extracted from the RefSeq of the 46 strains in GenBank. During similarity comparison and evolutionary tree analysis, multiple copies of 16S rRNA loci from one single mycoplasma strain were either retained with postfix, such as a, b, and c, to represent each copy or randomly selected and retained as only one copy. The 16S rRNA gene sequences and all single-copy core gene concatenated nucleotide sequences were aligned using MUSCLE integrated in MEGA X, respectively (26). The alignment results were used to construct phylogenetic trees using the maximum likelihood method based on the Tamura-Nei model in MEGA X with a bootstrap value of 1,000. Online software was used to calculate and draw a Venn diagram of the core-pan gene analysis results of six M. synoviae strains (http://bioinformatics.psb.ugent.be/webtools/Venn/). The homologies of complete genome sequences between every two M. synoviae strains were analyzed using BLASTN. The collinearity analysis and the number of SNPs (single nucleotide polymorphisms) and InDels (nucleotide insertion or deletion mutations) between two complete genome sequences of M. synoviae strains were analyzed using Mauve software (http://darlinglab.org/mauve/mauve.html). Insertion sequences (IS) were analyzed using the ISfinder tool (http://www-is.biotoul.fr).

The functional genes shared by the six strains of M. synoviae, and specific to the remaining 40 strains of Mycoplasma, were obtained from the gene set obtained by the core/pan gene clustering in the previous step. A BLASTP similarity search (https://blast.ncbi.nlm.nih.gov/Blast.cgi) was subsequently performed on these gene encoded amino acid sequences using default parameters against the non-redundant protein sequences (nr) database to obtain specific functional genes with no homology in the database except for M. synoviae. The functional genes that were ultimately obtained were used as candidates for molecular diagnostic target genes. The candidate target genes were subjected to a BLASTN similarity search (https://blast.ncbi.nlm.nih.gov/Blast.cgi) against the standard database. Gene coverage > 40% and an E value < 10⁻³ were used as cutoffs to identify genes specific to M. synoviae and absent in all other species in the database at the nucleotide level (11).

Primers targeting the genes of interest were designed using Primer3Plus (http://www.primer3plus.com/) and Primer Premier 5.00 (PREMIER Biosoft, http://www.premierbiosoft.com/). The specificity of the primers was determined using BLASTN against the standard databases. The qPCR assays were measured with the SYBR Premix Ex Taq II kit (TaKaRa, Japan) and run on the StepOnePlus^TM Real-Time PCR System (Applied Biosystems, US), according to the manufacturer's instructions.

Partial targeted gene sequences containing the qPCR primers were amplified by standard PCR using the primers designed as described above and cloned into the pMD-18T vector (TaKaRa, Japan) according to the user's manual. E. coli DH5α cells were used as the host for the recombinant plasmids. The copy number of recombinant plasmids was determined using an online calculator provided by the URI Genomics & Sequencing Center (http://cels.uri.edu/gsc/cndna.html) for determining the number of copies of a template. To evaluate sensitivity and reproducibility of the qPCR assays, the purified recombinant plasmids were diluted to 5 and 2.5 copies/reaction as well as continuously diluted by a factor of 10 to final concentrations of 10¹⁰ copies/reaction to 1 copy/reaction and used as the template to perform qPCR assays. Primers were also optimized by testing in the range of 100–500 nM.

To assess the specificity of the established qPCR method, the low pathogenic M. gallinarum (n = 1) and the major avian pathogenic bacteria/mycoplasmas, including M. iowae (n = 1), M. meleagridis (n = 1), M. synoviae (n = 20), M. gallisepticum (n = 3), S. aureus (n = 6), avian pathogenic E. coli (n = 4), S. Pullorum (n = 2), S. Gallinarum (n = 3), and S. Typhimurium (n = 3), were selected for analysis. The background information of these pathogens is provided in Supplementary Table 2. The last six pathogens mentioned above can cause arthritis in chickens. Clinically, it is often necessary to make a differential diagnosis of arthritis caused by these pathogens and M. synoviae (27–32). The genomic DNA of these 44 reference or clinical isolated strains of mycoplasma or bacteria was extracted. In addition, positive M. synoviae-infected chicken hock joint tissue samples, which were identified through culture-based pathogen isolation, genomic DNA extraction, and conventional PCR identification procedures reported in previous studies (33, 34), were selected for analysis. Chicken hock joint tissue samples positive for M. synoviae infection (n = 48) and negative control tissue samples (n = 10) were homogenized, and the mycoplasma genomic DNA was extracted to evaluate the specificity of the established qPCR method. Mycoplasma and bacterial genomic DNA was extracted using a Mycoplasma gDNA Mini Kit (BIOMIGA, China) and a MiniBEST Bacteria Genomic DNA Extraction Kit Ver.3.0 (TaKaRa, Japan), respectively, according to the manufacturers' instructions. All of the genomic and sample DNA concentrations were standardized to 1 ng/reaction. A sample was considered positive if the cycle threshold (C_T) value was <36 and negative if the C_T value was undetectable. Samples with C_T values equal to or >36 were re-prepared and the detection concentration was increased by 10 times to re-test.

The results of CD-HIT rapid clustering showed that there were 10,547 pan gene clusters among the 46 mycoplasma strains (Supplementary Table 3). After screening, a total of 16 protein-coding functional gene clusters (core gene clusters) were found to be shared by the 46 mycoplasma strains, of which three genes were multiple copies in some mycoplasmas and 13 genes were single copies in all 46 mycoplasmas (Table 1). Among these 16 genes, only gap was a double-copy in the six strains of M. synoviae, and the remaining 15 genes were single-copy in the six strains of M. synoviae (Table 1). In the pan gene clusters, 253 homologous gene clusters were shared by six M. synoviae strains and were absent in the other 40 mycoplasmas. Among these, 244 gene clusters were single-copy in all six strains of M. synoviae, and 9 gene clusters were multiple-copies in some of or all six M. synoviae strains (Supplementary Table 3).

We then performed genetic evolutionary analyses of M. synoviae relative to other major animal-derived mycoplasmas to study the specificity of M. synoviae. For homology comparison and evolutionary tree analysis, the nucleotide sequence of a gene is often more recognizable than the amino acid sequence it encodes. First, we used the DNA sequences of the 16S rRNA genes of each mycoplasma to perform homology comparisons and phylogenetic tree analysis. Given that most mycoplasmas have two or more copies of 16S rRNA, the phylogenetic tree based on the multi-copy 16S rRNA of all mycoplasma strains was not concise (Supplementary Figure 1). To facilitate the analysis, only one of the 16S rRNA multiple copies of each mycoplasma was randomly selected for phylogenetic tree analysis (Figure 1). The DNA sequences of 13 single-copy core genes from the 46 mycoplasma strains were subsequently selected for phylogenetic tree analysis (Figure 2). In general, the two phylogenetic trees were roughly the same.

Both phylogenetic tree analyses showed that M. synoviae was relatively distantly related to another major mycoplasma pathogen of chicken origin, M. gallisepticum. Phylogenetically, the two mycoplasmas were located in the two of the largest branches. The six strains of M. synoviae were placed into a relatively independent clade, and as a whole, were relatively closely related to a variety of mycoplasmas of avian origin, including Mycoplasma gallinaceum, Mycoplasma columborale, Mycoplasma glycophilum, Mycoplasma anatis, and Mycoplasma pullorum. Correspondingly, M. gallisepticum and M. iowae were closely related to two major human pathogenic mycoplasmas, Mycoplasma genitalium and Mycoplasma pneumoniae. In addition to M. gallisepticum and M. iowae, other avian mycoplasmas including M. synoviae, bovine mycoplasmas (except Mycoplasma mycoides), equine mycoplasma, and canine mycoplasmas were classified into one major clade.

The six strains of M. synoviae were divided into two and four category levels using 16S rRNA-based and 13 core gene-based phylogenetic tree analyses. The two evolutionary trees both showed that strains MS-H and 86079/7NS had higher homology between each other than the other four M. synoviae strains (Figures 1, 2). In the core gene-based phylogenetic tree analysis, the two subspecies of Mycoplasma capricolum were classified into one category, and then into the same category with M. mycoides. In the 16S rRNA-based phylogenetic tree analysis, the two subspecies of M. capricolum were not classified into the same lowest category. The three species of major swine pathogenic Mycoplasma were classified in the same category in the 13 core gene-based phylogenetic tree analyses, but not in the 16S rRNA-based analysis. In addition, M. meleagridis as well as canine mycoplasmas were classified into their own independent category in the 16S rRNA-based phylogenetic tree analysis, but not in the 13 core gene-based analysis.

In addition, the differences between M. synoviae and two species of pathogenic mycoplasmas that are characterized by causing mammalian arthritis and synovitis, Mycoplasma arthritidis and Mycoplasma hyosynoviae, were compared. Only two M. arthritidis strains (158L3-1 and NCTC10162) and one M. hyosynoviae strain (M60) had complete genome sequences and annotation information in GenBank. From the phylogenetic tree analysis, M. synoviae had low homology with M. arthritidis and M. hyosynoviae, while the homology relationship between M. arthritidis and M. hyosynoviae was very close. After screening the results of CD-HIT rapid clustering, there were 89 clusters in common between the six M. synoviae strains and two M. arthritidis strains, while none of the 89 clusters were absent in the rest of the 38 mycoplasma strains. There were 91 clusters in common between the six M. synoviae strains and one M. hyosynoviae strain, while none of the 91 clusters were absent in the rest of the 39 mycoplasma strains. There were 77 clusters in common among the six M. synoviae strains, two M. arthritidis strains, and one M. hyosynoviae strain, while none of the 77 clusters were absent in the rest of the 37 mycoplasma strains. Additionally, there was no cluster that was absent in MS-H but in common between five virulent M. synoviae strains and two M. arthritidis strains; between five virulent M. synoviae strains and one M. hyosynoviae strain; or among five M. synoviae strains, two M. arthritidis strains, and one M. hyosynoviae strain (Supplementary Table 3).

Collinearity analysis showed that the genomes of the six M. synoviae strains had basic collinearity. In addition to the previously reported MS-H and 86079/7NS, which have a unique inversion of approximately 55 kb relative to ATCC 25204 (10, 35), the same inversion was also found in NCTC10124 and HN01 (Figure 3A). The two ends of the inversion were bound by the coding genes of asparagine ligase (VY93_RS01590) and tRNA-Gly (VY93_RS01775). There was a P80 family lipoprotein (VY93_01545 and VY93_RS01790) on the two wings of the inversion.

Pairwise complete genome nucleotide sequence coverage, identity, SNP number, and InDel number of ATCC 25204 or MS-H with the other five M. synoviae strains ranged between 96.00 and 100.00%, 99.09 and 99.98%, 279 and 6,216, and 19 and 472, respectively (Table 2 and Supplementary Table 4). In general, the differences between ATCC 25204 and NCTC10124 and between MS-H and 86079/7NS were much smaller than between the other two strains listed in Table 2.

Core-pan gene analysis showed that there were 539 core clusters and 727 pan clusters in the six M. synoviae strains. Among them, only one cluster (cluster No. 10256, contained one gene, protein ID: WP_026365140.1, locus tag: MSH_RS02250, Supplementary Table 3) or two clusters that were unique to MS-H or 86079/7NS when compared to the other five strains, while 25 clusters were present in both MS-H and 86079/7NS but absent in the other four M. synoviae strains (Figure 4, Supplementary Table 5). These results indicate that MS-H and 86079/7NS have higher homology to each other relative to the other four strains. MSH_RS02250 encodes a hypothetical protein of 70 amino acids in MS-H (GenBank accession No. NZ_CP021129) and was subsequently subjected to a BLASTN similarity search. The result showed that the nucleotide sequence of MSH_RS02250 also exists in the corresponding position of the genomes of strains 53, ATCC 25204, 86079/7NS, and NCTC10124, whose GenBank accession Nos. are listed in Supplementary Table 1, but was not annotated as CDS in these strains. In addition, through SNP and InDel analysis, except in the vlhA locus, there were only point mutations and insertion/deletion mutations within 1 base in MS-H relative to its parent strain 86079/7NS (Supplementary Table 4). Therefore, MS-H had no specific functional genes, relative to the other five M. synoviae strains, that can be used for molecular diagnosis.

No cluster or seven clusters were unique to ATCC 25204 or NCTC10124 when compared to the other five strains, respectively, while 28 clusters were present in both ATCC 25204 and NCTC10124, but absent in the other four M. synoviae strains, indicating that ATCC 25204 and NCTC10124 have higher homology to each other relative to the other four strains. In contrast, strains 53 and HN01 had 18 and 22 unique clusters not found in the other five strains; while four clusters were present in both strains 53 and HN01 but absent in the other four M. synoviae strains, indicating that 53 and HN01 each have relatively distant homology with the other five strains (Figure 4).

A total of 253 clusters were found with common genes that were present in all six M. synoviae strains and absent in all of the other 40 mycoplasma strains, as mentioned above. Subsequently, 253 clusters containing genes were subjected to BLASTP similarity searches against the nr database in order to find unique functional genes in M. synoviae relative to species with known genome sequences. As a result, only six genes were found with amino acid sequences that could not be found by BLASTP in species other than M. synoviae (Table 3).

Two of the six M. synoviae-specific genes with the largest gene sizes, VY93_RS02250 and VY93_RS02300, are located between two housekeeping genes, Elongation factor Ts and P. The upstream and downstream regions of these two genes mostly consist of functional genes related to protein synthesis and energy supply, and no IS element was found by gene function annotation in GenBank or analysis using the ISfinder (Figure 3B and Supplementary Table 6). The sequences of VY93_RS02250, VY93_RS02300, and their upstream and downstream genes in the six strains of M. synoviae were consistent, and the DNA sequence homology of each gene was not <98.17% as found by BLASTN analysis (Supplementary Table 6). There was only one exception: the 53 strain had a single-base frameshift mutation at the position corresponding to the open reading frame (ORF) of ATCC 25204 VY93_RS02305, which caused the position to be divided into two genes, MS53_RS02170 and MS53_RS02175 (Supplementary Tables 4, 6). The DNA sequence homology of the two strains in the ORF where VY_RS02305 is located is 99.02%. BLASTN analysis showed that VY93_RS02250 and VY93_RS02300 were specific to M. synoviae, the coverage of the nucleotide sequences of these two genes in different M. synoviae strains was 100%, and the identity was >99% (Supplementary Table 6). These two genes were selected as targets for molecular diagnosis in this study. VY93_RS02250 and VY93_RS02300 were used as targets for establishing subsequent qPCR detection methods.

The primer sequences for the two qPCR methods are listed in Table 4. The product sizes of the qPCR assay detecting VY93_RS02250 and VY93_RS02300 were 102 and 173 bp, respectively. The optimal concentrations of all primers were determined to be 200 nM (data not shown). The efficiency of the two qPCR assays was tested on control plasmids with different DNA copy numbers (10¹⁰-10 copies/reaction). Drawing a standard curve with the lg numbers of DNA copy number as the abscissa and C_T values as the ordinate, the correlation between DNA copy numbers and C_T values for VY93_RS02250 was R² = 0.9997 (linear equation, y = −3.3393x + 38.267), and for VY93_RS02300 was R² = 0.9989 (linear equation, y = −3.3269x + 38.074).

The reproducibility of the two qPCR methods was evaluated with 10 replicates of each concentration of DNA copies in plasmids per reaction (Figures 5A,C). The results showed that when the number of DNA copies per reaction was 1,000, 100, and 10, the detected C_T values for genes VY93_RS02250 and VY93_RS02300 were 28.41 ± 0.38 and 28.62 ± 0.42, 31.78 ± 0.41 and 31.58 ± 0.32, and 34.75 ± 0.51 and 34.57 ± 0.52, respectively, which indicated good reproducibility of the two qPCR detection methods. Evaluation of the analytical sensitivity for the two targets showed 100% positive detection at a concentration of 10 DNA copies per reaction and variable detection percentage at a concentration of five or fewer DNA copies per reaction. Detection limits for VY93_RS02250 and VY93_RS02300 were 10 copies per reaction and 5 copies per reaction, respectively (Figures 5B,D).

The specificity test showed that when the genomic DNA of bacteria or mycoplasma strains were standardized to 1 ng/reaction, the positive amplified C_T values of the genes VY93_RS02250 and VY93_RS02300 of all 20 M. synoviae strains were 19.28 ± 0.95 and 19.60 ± 0.98, respectively; all other bacterial and mycoplasma strains had negative results (Supplementary Table 2). Clinical samples were also tested. When the primers targeting genes VY93_RS02250 and VY93_RS02300 were used for qPCR detection, and the extracted sample DNA standardized to 1 ng/reaction was used as the template, the test results showed that the C_T values of positive samples were ≤ 32.247, while the C_T value could not be detected in any of the negative samples (Supplementary Table 2).