# Pangenomic analysis of Coxiella burnetii unveils new traits in genome architecture

Rita Abou Abdallah, Matthieu Million, Jeremy Delerce, Hussein Anani, Awa Diop, Aurelia Caputo, Rita Zgheib, Elodie Rousset, Karim Sidi Boumedine, Didier Raoult, et al.

► To cite this version:

Rita Abou Abdallah, Matthieu Million, Jeremy Delerce, Hussein Anani, Awa Diop, et al. Pangenomic analysis of Coxiella burnetii unveils new traits in genome architecture. Frontiers in Microbiology, 2022, 13, pp.1022356. 10.3389/fmicb.2022.1022356. anses-03927880

HAL Id: anses-03927880

https://anses.hal.science/anses-03927880

Submitted on 6 Jan 2023

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
<figure><img src="image_2.png" /></figure>

Check for updates

OPEN ACCESS

EDITED BY

Vasco Ariston De Carvalho Azevedo,
Universidade Federal de Minas Gerais,
Brazil

REVIEWED BY

Sukhadeo B. Barbuddhe,
ICAR National Research Centre on
Meat (NRCM), India
Sandeep Tiwari,
Federal University of Minas Gerais,
Brazil

*CORRESPONDENCE

Pierre-Edouard Fournier
pierre-edouard.fournier@univ-amu.fr

SPECIALTY SECTION

This article was submitted to
Evolutionary and Genomic
Microbiology,
a section of the journal
Frontiers in Microbiology

RECEIVED 18 August 2022

ACCEPTED 10 October 2022

PUBLISHED 21 November 2022

CITATION

Abou Abdallah R, Million M, Delerce J,
Anani H, Diop A, Caputo A, Zgheib R,
Rousset E, Sidi Boumedine K, Raoult D
and Fournier P-E (2022) Pangenomic
analysis of *Coxiella burnetii* unveils
new traits in genome architecture.
*Front. Microbiol.* 13:1022356.
doi: 10.3389/fmicb.2022.1022356

COPYRIGHT

© 2022 Abou Abdallah, Million,
Delerce, Anani, Diop, Caputo, Zgheib,
Rousset, Sidi Boumedine, Raoult and
Fournier. This is an open-access article
distributed under the terms of the
Creative Commons Attribution License
(CC BY). The use, distribution or
reproduction in other forums is
permitted, provided the original
author(s) and the copyright owner(s)
are credited and that the original
publication in this journal is cited, in
accordance with accepted academic
practice. No use, distribution or
reproduction is permitted which does
not comply with these terms.

Pangenomic analysis of *Coxiella burnetii* unveils new traits in genome architecture

Rita Abou Abdallah<fnref n="1" />,², Matthieu Million²,³, Jeremy Delerce²,³,
Hussein Anani<fnref n="1" />,², Awa Diop<fnref n="1" />,², Aurelia Caputo²,³,
Rita Zgheib<fnref n="1" />,², Elodie Rousset⁴, Karim Sidi Boumedine⁴,
Didier Raoult²,³ and Pierre-Edouard Fournier<fnref n="1" />,²*

Coxiella burnetii is the etiological agent of Q fever, a worldwide zoonosis able to cause large outbreaks. The disease is polymorphic. Symptomatic primary infection is named acute Q fever and is associated with hepatitis, pneumonia, fever, and auto-immune complications while persistent focalized infections, mainly endocarditis, and vascular infections, occur in a minority of patients but are potentially lethal. In order to evaluate the genomic features, genetic diversity, evolution, as well as genetic determinants of antibiotic resistance, pathogenicity, and ability to cause outbreaks of Q fever, we performed a pangenomic analysis and genomic comparison of 75 C. burnetii strains including 63 newly sequenced genomes. Our analysis demonstrated that C. burnetii has an open pangenome, unique genes being found in many strains. In addition, pathogenicity islands were detected in all genomes. In consequence C. burnetii has a high genomic plasticity, higher than that of other intracellular bacteria. The core- and pan-genomes are made of 1,211 and 4,501 genes, respectively (ratio 0.27). The core gene-based phylogenetic analysis matched that obtained from multi-spacer typing and the distribution of plasmid types. Genomic characteristics were associated to clinical and epidemiological features. Some genotypes were associated to specific clinical forms and countries. MST1 genotype strains were associated to acute Q fever. A significant association was also found between clinical forms and plasmids. Strains harboring the QpRS plasmid were never found in acute Q fever and were only associated to persistent focalized infections. The QpDV and QpH1 plasmids were associated to acute Q fever. In addition, the Guyanese strain CB175, the most virulent strain to date, exhibited a unique MST genotype, a distinct COG profile and an important variation in gene number that may

<footnote n="1">¹Aix Marseille Université, Institut de Recherche pour le Développement (IRD), Service de Santé des Armées, AP-HM, UMR Vecteurs Infections Tropicales et Méditerranéennes (VITROME), Institut Hospitalo-Universitaire Méditerranée Infection, Marseille, France, ²Institut Hospitalo-Universitaire Méditerranée Infection, Marseille, France, ³Aix-Marseille Université, Institut de Recherche pour le Développement (IRD), UMR Microbes Evolution Phylogeny and Infections (MEPHI), Institut Hospitalo-Universitaire Méditerranée Infection, Marseille, France, ⁴French Agency for Food, Environmental and Occupational Health Safety (ANSES), Sophia Antipolis Laboratory, Animal Q Fever Unit, Sophia Antipolis, France</footnote>
explain its unique pathogenesis. Therefore, strain-specific factors play an important role in determining the epidemiological and clinical manifestations of Q fever alongside with host-specific factors (valvular and vascular defects notably).

KEYWORDS

Coxiella burnetii, Q fever, pangenome, genotyping, pathogenesis

Introduction

Coxiella burnetii, an obligate Gram-negative and strictly
intracellular bacterium responsible for the zoonotic disease Q
fever (Angelakis and Raoult, 2010), was first discovered in 1935
during an outbreak of febrile illness in abattoir workers in
Brisbane, Australia (Derrick, 1983). Two years later, C. burnetii
was first isolated after injecting patient blood and urine samples
to guinea pigs (Angelakis and Raoult, 2010). It was initially
classified as a member of the genus Rickettsia (Angelakis and
Raoult, 2010), prior to being reclassified into the new genus
Coxiella in 1948, after studying its cultural and biochemical
characteristics. Currently, Q fever is endemic worldwide, except
in New Zealand (Hilbink et al., 1993). The disease occurs
mostly as sporadic infections, although outbreaks have been
reported in France (Frankel et al., 2011), Great Britain (Smith
et al., 1993), Italy (Manfredi Selvaggi et al., 1996), Germany
(Lyytikäinen et al., 1998), and the Netherlands (de Bruin et al.,
2012). The largest outbreak occurred from 2007 to 2012 in the
latter country, with more than 4,000 diagnosed cases (de Bruin
et al., 2012).

Even though its animal reservoir is large, domestic
ruminants such as goats, sheep and less frequently cattle,
represent the most common source of Q fever. *Coxiella burnetii*
is mostly transmitted to humans by aerosols from parturition
fluids of infected animals or from microorganisms spread by the
wind. The bacterium is very resistant and may survive several
weeks in the environment. Although *C. burnetii* was detected
in ticks, this mode of transmission remains rare (*Angelakis and
Raoult, 2010*).

In animals, Q fever is almost asymptomatic with the
exception of cases of abortion, low birth weight and infertility
(Angelakis and Raoult, 2010). In contrast to animal infection,
human Q fever is characterized by a polymorphic clinical
presentation (Raoult et al., 2005; Eldin et al., 2017; Melenotte
et al., 2020). In acute, or primary Q fever, symptoms can
range from an asymptomatic sero-conversion to a flu-like self-
limited illness, to more severe symptoms including interstitial
pneumonia (Caron et al., 1998), granulomatous hepatitis and
rarely cardiac involvement (Fournier et al., 2001), skin rash
and neurological signs (Angelakis et al., 2012). Patients may
subsequently develop persistent focalized infections defined by

both the microbe multiplication and the presence of organic
lesions (Raoult, 2017). Endocarditis is the most frequent
persistent focalized infection with C. burnetii (Angelakis and
Raoult, 2010) involving in priority the aortic and mitral valves
as well as valvular prostheses. Vascular infections are the second
most frequent persistent focalized infection (O'Donnell et al.,
2007; Angelakis and Raoult, 2010; Million and Raoult, 2017)
and have a more severe prognosis than endocarditis. Other, less
frequent manifestations of persistent focalized Q fever include
osteoarticular infections (Angelakis et al., 2014), pseudo-tumors
of the lung and persistent lymphadenitis that may result in
lymphoma (Melenotte et al., 2016). It was proposed that the
clinical presentation, host type and severity of Q fever may vary
depending on the involved C. burnetii strain (Glazunova et al.,
2005; Roest et al., 2011; de Bruin et al., 2012; Angelakis et al.,
2013).

The availability of whole C. burnetii genome sequences
was a major step forward for studying the pathophysiology of
Q fever (Glazunova et al., 2005; D'Amato et al., 2014, 2016;
Kuley et al., 2017). The first C. burnetii genome was sequenced
from Nine Mile Phase I RSA493 (NM-I) strain in 2003. The
1,995,275-bp-long NM-I genome with a 37,393 QpH1 plasmid
had mobile elements and had different metabolic capabilities,
transport profile and had limited extent of genome reduction
compared to obligate intracellular pathogens like Rickettsia,
Chlamydia, and Mycobacterium leprae suggesting comparatively
recent intracellular adaptation (Seshadri et al., 2003). Over
recent years, thanks to the development of high throughput
sequencing methods, 14 additional C. burnetii strain genomes
were sequenced and allowed comparative genomic analyses.
Prior to our study, the size and G + C content of C. burnetii
genomes ranged from 1,956,650 to 2,093,477 bp and from
42.5 to 42.9%, respectively (Kuley et al., 2017). Four plasmid
types were also identified (QpH1, QpRs, QpDv, and QpDg),
with strains exhibiting any one or none of them (Beare et al.,
2009; Sidi-Boumedine et al., 2014). Samuel et al. (1985) studied
the correlation between plasmid type and disease presentation
(acute or persistent focalized) and suggested that plasmid types
may be associated to distinct clinical forms. Plasmid type has
also been associated with pregnancy outcome (Angelakis et al.,
2013). Beare et al. (2009) proposed that C. burnetii strains
are at various stages of pathoadaptation. They characterized
the genomic features associated to the obligate intracellular lifestyle of this pathogen such as genome reduction, presence of pseudogenes, and mobile genetic elements (Beare et al., 2009). The role of up-regulated genes in the survival and replication of the bacterium in hosts was also investigated (Kuley et al., 2015). These authors suggested that the Dugway strain was the most primitive due to its large genome and few pseudogenes and insertion sequences, in addition to the fact that it has the shortest distance to the root of the phylogenetic tree (Beare et al., 2009). By studying the virulent strain Guyane CB175, D'Amato et al. (2015) showed that its genome exhibited a 6,105 bp deletion and a large number of gene mutations compared to strain NM-1. Other genomic comparisons aimed at deciphering the high virulence of some strains such as Guyane CB175 and those responsible for the German outbreaks (D'Amato et al., 2014, 2015; Kuley et al., 2017). As whole-genome sequencing and comparative genome analysis of a large number of related strains has recently emerged as a cost-effective and convenient method allowing a better identification of the diversity and composition of the global gene repertoire of a bacterium (Tettelin et al., 2008), we performed a pangenomic analysis of C. burnetii using the genomes from 75 C. burnetii strains including 63 novel genomes, in order to evaluate the genetic diversity, evolution, antibiotic resistance, and disease-specific characteristics of this pathogen.

# Materials and methods

## Bacterial strains

A total of 75 strains were included in this retrospective study. Sixty-three strains had been isolated from clinical specimens of patients suffering from acute or persistent focalized Q fever. No patient objected the use of their strain for research purpose. Another 12 strains were isolated from ticks or animals. Among the 75 studied genomes, eight were obtained from Genbank (strains RSA493, RSA331, Dugway, CbuK_Q154, CbuG_Q212, 3262, MSUGoat_Q177, Z3055). All other genomes had been sequenced and assembled in our laboratory. Among those sequenced in our laboratory, the genomes of four strains were previously published (CB175, CB196, CB185, Schperling).

## Identification of *Coxiella burnetii* strains, growth conditions, and DNA extraction

*Coxiella burnetii* strains isolated in our laboratory were identified using PCR targeting the IS1111 and IS30 repeated insertion sequences. All strains were cultured, grown, and produced on Vero cell monolayers (ATCC CRL 1587) maintained in Minimal Essential Medium (MEM) (Invitrogen,

Cergy-Pontoise, France) supplemented with 4% fetal bovine serum (Invitrogen) and 1% L-glutamine (Invitrogen), in a 5% CO₂ atmosphere at 35°C. After cytopathic effect detection, bacteria were purified by successive centrifugations in order to remove cells. DNA extraction was performed using the automated extraction EZ1 DNA Tissue Kit (Qiagen, Hilden, Germany) according to the manufacturer's instructions.

## Whole genome sequencing

To investigate the clinical features of C. burnetii, we sequenced the genomes of human isolates causing various clinical symptoms, in addition to five strains isolated from ticks and six from animals. Among 67 strains sequenced in our laboratory, 15 had been sequenced using SOLiD 4_Life technology while the remaining 52 were sequenced using the Illumina Miseq technology. For the SOLiD 4_Life sequencing, the paired-end library was constructed from 1 µg of purified genomic DNA. The DNA was fragmented on the Covaris device. Libraries were pooled in equimolar ratios and size selected on the E-Gel iBase system at 240–270 bp. Sequencing was carried out to 50 × 35 bp using SOLiD™ V4 chemistry on one full slide of 96 barcoded projects on an Applied Biosystems SOLiD 4 sequencer. The output read length was as expected 85 bp.

The paired-end strategy was used for the 52 strains sequenced using MiSeq Technology (Illumina Inc., San Diego, CA, USA). The genomic DNA was quantified by a Qubit assay with the high sensitivity kit (Life Technologies, Carlsbad, CA, USA) and was barcoded in order to be mixed with 14 other projects with the Nextera XT DNA sample prep kit (Illumina Inc., San Diego, CA, USA).

Purified, normalized and barcoded DNA were pooled for sequencing on the MiSeq sequencer. The pooled single strand library was loaded onto the reagent cartridge and then onto the instrument along with the flow cell. Automated cluster generation and paired-end sequencing with dual index reads were performed in a single 39-h run in 2 × 250-bp.

## Genome assembly and functional annotation

Prior to assembly, raw data quality was checked with the fastQC tool version 0.11.8¹ (Supplementary Table 1).

Genomes sequenced with Solid were assembled by mapping using the CLC Genomics Workbench version 7.5 (Qiagen). Genomes sequenced with Illumina Miseq were assembled using SPAdes version 3.10 (Bankevich et al., 2012). All C. burnetii genomes were annotated using Prokka version 1.12-beta (Seemann, 2014).

<footnote n="1">1 https://www.bioinformatics.babraham.ac.uk/projects/fastqc/</footnote>
## Genotyping of *Coxiella burnetii* strains

The genotypes from all studied strains were determined using the multispacer sequence typing (MST) method, based on 10 intergenic regions as described previously (Glazunova et al., 2005). Moreover, MST PCR was carried out for new or doubtful MST genotypes (Glazunova et al., 2005).

## Comparison of *Coxiella burnetii* genome contents

We performed a COG (Cluster of Orthologous Groups) analysis on the 75 whole genomes as well as on the core genome. We looked for insertion sequences using the ISfinder web server. We also explored the presence of resistance genes using the RAST server. Moreover, digital DNA-DNA hybridization (dDDH) and OrthoAni analysis were performed on the 75 genomes in order to evaluate the range of genomic similarity between strains (Auch et al., 2010; Lee et al., 2016). Plasmids were determined in silico with blast and mapping. Virulence factors were determined using the VFDB database (Chen, 2004). The IslandViewer web server was used to predict genomic islands (Bertelli et al., 2017).

## Phylogenomic and pangenomic analysis

Pangenome analyses along with core and non-core region analyses were performed with Proteinortho version 5, with an amino acid identity percentage and coverage threshold of 40 and 60%, respectively (Lechner et al., 2011). Prokka annotations were provided to Proteinortho as an input. Briefly, core regions were defined as genome sequences present in all 75 genomes, while the non-core regions were defined as those differentially present in studied genomes. Unique proteins in sequenced strains were identified from the orthologous list generated by Proteinortho after manual correction of split genes. We constructed a phylogenomic tree using the core region of the 75 C. burnetii genomes rather than few genomic loci. Core-genome CDS were inferred from the pangenome, concatenated and aligned using a custom Python script, and then used to construct a phylogenetic tree with MEGA6. SNPs were extracted from the core genome using the SNP-sites software (Page et al., 2016). Visualizations of phylogenetic trees were produced with FigTree.

## Statistical analysis

Statistical analysis was performed using the SPSS software version 22 and XLSTAT 2019 (Addinsoft, Paris, France).

## Correlation between genomic features of the 75 genomes

The Spearman correlation coefficient was used to correlate between quantitative variables. We studied the relations between the presence of an efflux pump, the genome size and G + C content.

## Quantitative principal component analysis

In order to explore the geographic and clinical data that corresponded to the 75 genomes, we conducted a principal component analysis with genomic information (genome size, G + C content, MST and plasmid type), geography, host, and clinical presentation.

## Logistic (exponential) principal component analysis

Thanks to the 75 genomes sequenced here and the data already analyzed in one of our previous studies (Glazunova et al., 2005; Angelakis et al., 2013), we completed the study (205 strains, Supplementary Table 2) of the association of genotypes defined by MST or plasmids and the country of origin of the strain or clinical presentation (acute Q fever or persistent focalized infection). Excluding quantitative variables (genome size, G + C content) and including only strains without missing data into a logistic (exponential) principal component analysis specifically adapted to the qualitative data (Figure 6).

# Results

## Genome sequences of *Coxiella burnetii* strains

In the present study, we reported the new genome sequences from 63 C. burnetii strains.

## Genomic features

The average genome size of C. burnetii isolates was 2.006 Mb (ranging from 1.910 to 2.158 Mb), and the mean G + C content was 42.62% (ranging from 42.4 to 42.9%), with a number of contigs ranging from 1 to 282 (Table 1).

## Genotyping

The 75 studied isolates belonged to 22 different MST genotypes. One new MST genotype was identified in this study (MST 60), and the first genome of a pre-identified MST was sequenced (MST55). The most frequent MST genotypes were MST 8 and 33 (11 and 8 strains, respectively). In contrast, the less represented genotypes were MST 19, 3, 23, 55, and 60, with only one strain each.
<table><caption>TABLE 1 General features of the 75 studied *Coxiella burnetii* genomes.</caption><thead><tr><th>Strain</th><th>G + C content</th><th>Genome size (pb)</th><th>No. of contigs</th><th>Accession number</th><th>Country</th><th>Source</th><th>Acute/ Persistent</th><th>Disease</th><th>MST genotype</th></tr></thead><tbody><tr><td>3262</td><td>42.9</td><td>2,093,477</td><td>1</td><td>CP013667</td><td>The Netherlands</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>BRASOV</td><td>42.7</td><td>1,910,413</td><td>70</td><td>CP103435</td><td>Romania</td><td>Human</td><td>Acute</td><td>Pneumonia</td><td>16</td></tr><tr><td>CB181</td><td>42.6</td><td>1,972,104</td><td>42</td><td>JANTNR000000000</td><td>French Guiana</td><td>Human</td><td>Acute</td><td>Pneumonia</td><td>17</td></tr><tr><td>CB1</td><td>42.6</td><td>2,009,715</td><td>78</td><td>JAOXFE000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>1</td></tr><tr><td>CB10</td><td>42.5</td><td>2,041,860</td><td>86</td><td>JAOXFD000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>8</td></tr><tr><td>CB111</td><td>42.6</td><td>2,019,099</td><td>72</td><td>JAOXFC000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>6</td></tr><tr><td>CB118</td><td>42.6</td><td>2,002,150</td><td>75</td><td>JAOXFB000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>3</td></tr><tr><td>CB119</td><td>42.6</td><td>1,981,226</td><td>47</td><td>JAOXFA000000000</td><td>Senegal</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>19</td></tr><tr><td>CB121</td><td>42.6</td><td>1,923,427</td><td>157</td><td>CP103434</td><td>Switzerland</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>12</td></tr><tr><td>CB125</td><td>42.5</td><td>2,040,530</td><td>76</td><td>JAOXEZ000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>8</td></tr><tr><td>CB13</td><td>42.6</td><td>1,933,670</td><td>142</td><td>CP103433</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>11</td></tr><tr><td>CB132</td><td>42.6</td><td>1,934,159</td><td>130</td><td>CP103432</td><td>Belgium</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>11</td></tr><tr><td>CB144</td><td>42.6</td><td>1,956,889</td><td>120</td><td>JAOXEY000000000</td><td>Senegal</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>12</td></tr><tr><td>CB149</td><td>42.7</td><td>1,924,741</td><td>152</td><td>CP103431</td><td>Senegal</td><td>Tick</td><td>NA</td><td>NA</td><td>6</td></tr><tr><td>CB152</td><td>42.5</td><td>2,030,144</td><td>51</td><td>JAOXEX000000000</td><td>Senegal</td><td>Tick</td><td>NA</td><td>NA</td><td>16</td></tr><tr><td>CB155</td><td>42.6</td><td>1,924,716</td><td>155</td><td>CP103430</td><td>United Kingdom</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>55</td></tr><tr><td>CB156</td><td>42.6</td><td>2,020,453</td><td>71</td><td>JAOXEW000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>4</td></tr><tr><td>CB157</td><td>42.8</td><td>2,022,809</td><td>113</td><td>JAOXEV000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB162</td><td>42.5</td><td>2,017,686</td><td>72</td><td>JAOXEU000000000</td><td>USA</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>8</td></tr><tr><td>CB170</td><td>42.6</td><td>1,929,924</td><td>151</td><td>CP103429</td><td>USA</td><td>Human</td><td>Persistent</td><td>Bone infection</td><td>8</td></tr><tr><td>CB171</td><td>42.6</td><td>2,050,297</td><td>94</td><td>JAOXET000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>4</td></tr><tr><td>CB175</td><td>42.6</td><td>1,989,609</td><td>1</td><td>HG825990</td><td>French Guiana</td><td>Human</td><td>Others</td><td>Endocarditis</td><td>17</td></tr><tr><td>CB176</td><td>42.6</td><td>1,970,312</td><td>37</td><td>JAOXES000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>17</td></tr><tr><td>CB180</td><td>42.6</td><td>1,944,582</td><td>112</td><td>CP103428</td><td>France</td><td>Human</td><td>Persistent</td><td>Bone infection</td><td>18</td></tr><tr><td>CB184</td><td>42.8</td><td>2,026,895</td><td>120</td><td>JAOXER000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>32</td></tr><tr><td>CB185</td><td>42.6</td><td>1,989,218</td><td>282</td><td>GCA_000470495</td><td>France</td><td>Human</td><td>Persistent</td><td>Miscarriage</td><td>33</td></tr><tr><td>CB190</td><td>42.5</td><td>2,018,107</td><td>68</td><td>JAOXEQ000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>8</td></tr><tr><td>CB192</td><td>42.5</td><td>2,020,838</td><td>69</td><td>JAOXEP000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>8</td></tr><tr><td>CB193</td><td>42.6</td><td>2,025,795</td><td>69</td><td>JAOXEO000000000</td><td>French Guiana</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>8</td></tr><tr><td>CB196</td><td>42.6</td><td>2,009,497</td><td>72</td><td>GCA_000820465</td><td>Saudi Arabia</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>4</td></tr><tr><td>CB198</td><td>42.6</td><td>2,003,740</td><td>77</td><td>JAOXEN000000000</td><td>France</td><td>Human</td><td>Acute</td><td>Pericarditis</td><td>1</td></tr><tr><td>CB2</td><td>42.8</td><td>2,031,973</td><td>117</td><td>JAOXEM000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>13</td></tr><tr><td>CB201</td><td>42.6</td><td>2,021,166</td><td>70</td><td>JAOXEL000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>4</td></tr><tr><td>CB202</td><td>42.5</td><td>2,009,975</td><td>86</td><td>JAOXEK000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>6</td></tr><tr><td>CB204</td><td>42.6</td><td>2,001,464</td><td>76</td><td>JAOXEJ000000000</td><td>France</td><td>Human</td><td>Acute</td><td>Isolated fever</td><td>1</td></tr><tr><td>CB207</td><td>42.4</td><td>2,118,168</td><td>274</td><td>JAOXEI000000000</td><td>France</td><td>Human</td><td>Acute</td><td>Hepatitis</td><td>16</td></tr><tr><td>CB211</td><td>42.5</td><td>2,046,524</td><td>55</td><td>JAOXEH000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>16</td></tr><tr><td>CB212A</td><td>42.7</td><td>1,995,348</td><td>71</td><td>JAOXEG000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>16</td></tr><tr><td>CB212</td><td>42.6</td><td>1,967,025</td><td>40</td><td>JAOXEF000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Aneurysm infection</td><td>20</td></tr><tr><td>CB213</td><td>42.8</td><td>2,040,940</td><td>115</td><td>JAOXEE000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>32</td></tr><tr><td>CB215</td><td>42.9</td><td>2,058,005</td><td>114</td><td>JAOXED000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>13</td></tr><tr><td>CB216</td><td>42.8</td><td>2,031,442</td><td>114</td><td>JAOXEC000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB217</td><td>42.9</td><td>2,060,311</td><td>118</td><td>JAOXEB000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB218</td><td>42.6</td><td>1,968,958</td><td>38</td><td>JAOXEA000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>20</td></tr><tr><td>CB220</td><td>42.6</td><td>1,969,527</td><td>40</td><td>JAOXDZ000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Bone infection</td><td>20</td></tr><tr><td>CB221</td><td>42.6</td><td>1,969,946</td><td>32</td><td>JAOXDY000000000</td><td>France</td><td>Tick</td><td>NA</td><td>NA</td><td>16</td></tr><tr><td>CB222</td><td>42.5</td><td>2,018,652</td><td>91</td><td>JAOXDX000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>6</td></tr><tr><td>CB3</td><td>42.6</td><td>2,002,831</td><td>77</td><td>JAOXDW000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>1</td></tr></tbody></table>

(Continued)
<table><caption>TABLE 1 (Continued)</caption><thead><tr><th>Strain</th><th>G + C content</th><th>Genome size (pb)</th><th>No. of contigs</th><th>Accession number</th><th>Country</th><th>Source</th><th>Acute/ Persistent</th><th>Disease</th><th>MST genotype</th></tr></thead><tbody><tr><td>CB4</td><td>42.5</td><td>2,016,251</td><td>60</td><td>JAOXDV000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>21</td></tr><tr><td>CB5</td><td>42.8</td><td>2,023,365</td><td>115</td><td>JAOXDU000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>12</td></tr><tr><td>CB6</td><td>42.9</td><td>2,050,803</td><td>129</td><td>JAOXDT000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>12</td></tr><tr><td>CB65</td><td>42.6</td><td>2,038,254</td><td>80</td><td>JAOXDS000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>9</td></tr><tr><td>CB68</td><td>42.6</td><td>1,928,614</td><td>149</td><td>CP103427</td><td>France</td><td>Mouse</td><td>NA</td><td>NA</td><td>9</td></tr><tr><td>CB69</td><td>42.9</td><td>2,051,941</td><td>115</td><td>JAOXDR000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>13</td></tr><tr><td>CB74</td><td>42.7</td><td>2,004,352</td><td>114</td><td>JAOXDQ000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB80</td><td>42.6</td><td>1,928,715</td><td>192</td><td>CP103426</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB85</td><td>42.6</td><td>1,957,952</td><td>118</td><td>JAOXDP000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>33</td></tr><tr><td>CB9</td><td>42.5</td><td>2,024,753</td><td>74</td><td>JAOXDO000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>8</td></tr><tr><td>CB94</td><td>42.4</td><td>2,068,894</td><td>176</td><td>JAOXDN000000000</td><td>France</td><td>Human</td><td>Acute</td><td>Hepatitis</td><td>16</td></tr><tr><td>CbuG-Q212</td><td>42.6</td><td>2,008,870</td><td>1</td><td>CP001019</td><td>Canada</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>21</td></tr><tr><td>CbuK-Q154</td><td>42.7</td><td>2,063,100</td><td>1</td><td>CP001020</td><td>USA</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>8</td></tr><tr><td>DogUtad</td><td>42.5</td><td>1,956,891</td><td>76</td><td>CP103425</td><td>Canada</td><td>Animal</td><td>Dog</td><td>NA</td><td>21</td></tr><tr><td>Dugway</td><td>42.4</td><td>2,158,758</td><td>1</td><td>CP000733</td><td>USA</td><td>Animal</td><td>Rodent</td><td>NA</td><td>60</td></tr><tr><td>IXO</td><td>42.6</td><td>1,974,499</td><td>67</td><td>JAOXDM000000000</td><td>Slovakia</td><td>Tick</td><td>NA</td><td>NA</td><td>23</td></tr><tr><td>CB225</td><td>42.6</td><td>2,002,796</td><td>78</td><td>JAOXDL000000000</td><td>France</td><td>Human</td><td>Acute</td><td>Hepatitis</td><td>1</td></tr><tr><td>CB220B</td><td>42.6</td><td>1,972,402</td><td>47</td><td>JAOXDK000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Bone infection</td><td>20</td></tr><tr><td>MSU_Q177</td><td>42.7</td><td>2,063,021</td><td>1</td><td>CP018150</td><td>USA</td><td>Goat</td><td>NA</td><td>NA</td><td>8</td></tr><tr><td>NM-I</td><td>42.7</td><td>1,995,281</td><td>1</td><td>NC_002971</td><td>USA</td><td>Tick</td><td>NA</td><td>NA</td><td>16</td></tr><tr><td>CB224</td><td>42.6</td><td>2,032,061</td><td>72</td><td>JAOXDJ000000000</td><td>France</td><td>Human</td><td>Persistent</td><td>Bone infection</td><td>8</td></tr><tr><td>PokerCat</td><td>42.5</td><td>2,030,687</td><td>68</td><td>JAOXDI000000000</td><td>Canada</td><td>Cat</td><td>NA</td><td>NA</td><td>21</td></tr><tr><td>Q229</td><td>42.5</td><td>1,975,920</td><td>41</td><td>JAOXDH000000000</td><td>Canada</td><td>Human</td><td>Persistent</td><td>Endocarditis</td><td>21</td></tr><tr><td>Henzerling RSA331</td><td>42.8</td><td>2,016,427</td><td>1</td><td>CP014559</td><td>Italy</td><td>Human</td><td>Acute</td><td>Pneumonia</td><td>18</td></tr><tr><td>Schperling</td><td>42.6</td><td>2,004,282</td><td>73</td><td>GCA_002634065</td><td>Kyrgyzstan</td><td>Human</td><td>Chronic</td><td>Fever</td><td>2</td></tr><tr><td>Z2534</td><td>42.9</td><td>2,056,819</td><td>115</td><td>JAOXDG000000000</td><td>Austria</td><td>Goat</td><td>NA</td><td>NA</td><td>32</td></tr><tr><td>Z3055</td><td>42.6</td><td>1,995,457</td><td>1</td><td>LK937696</td><td>Germany</td><td>Ewe</td><td>NA</td><td>NA</td><td>33</td></tr></tbody></table>

## Plasmid analysis and resistance genes of *Coxiella burnetii* isolates

Four types of plasmids were found in 70 out of the 75 studied genomes, while five strains were plasmidless (Figure 1). Among these 70 strains, 40, 18, and 11 strains exhibited the QpH1, QpRS, and QpDv plasmids, respectively. In addition, the Dugway strain was the only strain to exhibit the QpDg plasmid.

Beta-lactam and fluoroquinolone resistance-coding genes were found in all strains. Nineteen strains exhibited a multidrug efflux pump coupled to Na⁺ and belonging to the MATE family. Regarding virulence factors, genes coding for the secretion system type IVB, were found in the 75 studied genomes.

## Phylogenomic analysis

Based on the clustering analysis, the strains could be divided into 13 clusters (Figure 1). The number of strains detected in each phylogenetic cluster along with country of detection,

MST genotype, exhibited plasmid and type of infection is summarized in Table 2. This phylogenetic tree is consistent with that obtained from MST sequence analysis (Supplementary Figure 1). Strains exhibiting the same plasmid type were grouped together. In addition, plasmidless strains were clustered into one isolated group (MST21).

## Comparative genome analysis and genomic characteristics

Comparative genome analyses were performed for all genomes. In order to determine the range of intra-species genomic similarity, we used dDDH and OrthoANI among strain pairs. dDDH values ranged from 74.5% between strains CB85 and CB3262, to 100% (Supplementary Table 3). OrthoANI values ranged from 99.37% between strain Z3055 and CB65, to 100% (Supplementary Table 4). The 100% values were always observed between strains from the same MST genotype and plasmid type.
<figure><img src="image_3.png" /><figcaption>FIGURE 1</figcaption></figure>

Phylogenetic tree of 75 *Coxiella burnetii* strains using core genome regions, using the Maximum Likelihood method within MEGA6. Branch color indicates the associated disease. Blue branch: persistent focalized Q fever, yellow branch: acute Q fever, green branch: tick isolate, purple branch: animal isolate, black branch: others, gray branch: no data available. *: presence of an efflux pump. On the right panel: the first column represents the MST genotype, the remaining columns represent the presence of plasmids with each symbol indicating a plasmid type.
<table><caption>TABLE 2 Summary of phylogenetic analysis of 75 *Coxiella* genomes.</caption><thead><tr><th>Phylogenetic cluster</th><th>No. of strains</th><th>Countries</th><th>MST genotypes</th><th>Plasmid</th><th>Infection types</th><th>Clinical manifestation</th></tr></thead><tbody><tr><td>1</td><td>20</td><td>France, Belgium, Switzerland, Austria, Germany, and Senegal</td><td>11, 12, 13, 32, 33</td><td>QpH1</td><td>Persistent</td><td>Primarily endocarditis</td></tr><tr><td>2</td><td>4</td><td>France, Italy, and Slovakia</td><td>18, 23</td><td>QpH1</td><td>Persistent</td><td>Pneumonia and bone infection</td></tr><tr><td>3</td><td>4</td><td>France</td><td>20</td><td>QpH1</td><td>Persistent</td><td>Bone infection, aneurysm, endocarditis</td></tr><tr><td>4</td><td>1</td><td>Senegal</td><td>19</td><td>QpH1</td><td>Persistent</td><td>Endocarditis</td></tr><tr><td>5</td><td>10</td><td>France, French Guiana, Romania, USA, and Senegal</td><td>16, 17</td><td>QpH1</td><td>Acute and Persistent</td><td>Endocarditis, hepatitis, pneumonia, aneurysm</td></tr><tr><td>6</td><td>1</td><td>France</td><td>16</td><td>QpH1</td><td>Acute</td><td>Hepatitis</td></tr><tr><td>7</td><td>5</td><td>Canada France</td><td>21</td><td>None</td><td>Persistent</td><td>Endocarditis</td></tr><tr><td>8</td><td>11</td><td>France, French Guiana, and the USA</td><td>8</td><td>QpRs</td><td>Persistent</td><td>Bone infection, aneurysm, endocarditis</td></tr><tr><td>9</td><td>2</td><td>France</td><td>9</td><td>QpRs</td><td>Persistent</td><td>Endocarditis</td></tr><tr><td>10</td><td>1</td><td>UK</td><td>55</td><td>QpRs</td><td>Persistent</td><td>Endocarditis</td></tr><tr><td>11</td><td>4</td><td>France, Senegal</td><td>6</td><td>QpRs</td><td>Persistent</td><td>Endocarditis</td></tr><tr><td>12</td><td>11</td><td>France, Saudi Arabia, and Italy</td><td>1, 2, 3, 4</td><td>QpDv</td><td>Acute and persistent</td><td>Endocarditis, hepatitis, isolated fever, pericarditis, fever, aneurysm infection</td></tr><tr><td>13</td><td>1</td><td>USA (Rodent)</td><td>60</td><td>QpDg</td><td>NA</td><td>NA</td></tr></tbody></table>

Regarding the insertion sequences, the 75 genomes contained the IS1111C, IS1111A, and ISCbu1 insertion sequences. As most genomes were draft sequences, we were not able to determine the number of IS sequences.

The COG analysis showed a homogenous distribution of gene categories among most genomes (Figure 2). Seventeen genomes had significantly higher numbers of genes in three COG categories. Strain CB175 (strain Guyane) showed the highest number of genes classified in the C to Z COGs categories (Figure 2). Strains CB207 and CB94 had a similar COG profiles with higher numbers of genes in the G to Z COG categories. Similarly, strains CB213, Z2534, CB217, CB215, CB69, CB6, 3262, CB2, CB216, CB184, Dugway, CB5, CB157, and CB74 had higher numbers of genes in L to Z categories (Figure 2).

## *Coxiella burnetii* pangenome

The pangenomic analysis identified a total of 4,501 genes, including 1,211 core genes, and 3,290 accessory genes, including 1,344 unique genes that were only present in one strain. The number of unique genes per strain ranged from 7 to 212. The CB175 and Dugway strains had the highest numbers of unique genes (212 and 140, respectively). Among unique genes, five were likely to be acquired from members of the phylum *Chlamydiae*. Strain CB5 presented three unique genes found in *Chlamydia abortus*, a species that causes fetal death in mammals. Of these genes, one coded a transketolase, and the other two are hypothetical proteins. Strain CB193 exhibited a unique gene coding an hypothetical protein also present in *Chlamydia*

*trachomatis*. Strain Dugway presented one gene found in *Parachlamydia acanthamoebae* and coding an aminoglycoside phosphotransferase, an enzyme known to cause resistance to this antibiotic family.

The core genome length was 1,789,492 bp with a G + C content of 42.7% and 10,192 SNP's. The variations of core genes, core genome/pangenome ratio and of the total numbers of genes with the addition of new genome sequences is represented in Figure 3. Within each of the above-defined phylogenetic clusters (Figure 1), the contribution of the core genome to the complete gene repertoire of sequenced *C. burnetii* strains ranged from 72 to 95%. However, when this analysis was performed on the 75 genomes, this core genome/pangenome percentage dropped to 27%. Figure 3 shows that almost all the important core/pangenome percentage shifts occurred between MST groups, in addition to an important variation for CB175. The presence/absence gene matrix (Figure 4) showed similar profiles in strains belonging to a given MST genotype. Genome plasticity was clearly between different MST genotypes. Genomic islands which are defined as clusters of genes of probable horizontal origin were found in all strains and in the virulent CB175 Guyanese strain.

## Statistical analysis

### Correlation between genomic features of the 75 genomes

A significant correlation was found between the presence of the MATE family efflux pump and the G + C content
<figure><img src="image_4.png" /><figcaption>FIGURE 2</figcaption></figure>

Distribution of genes into cluster of orthologous groups (COG) categories. J: translation, A: RNA processing and modification, K: transcription, L: Replication, recombination and repair, B: Chromatin structure and dynamics, D: Cell cycle control, mitosis and meiosis, Y: Nuclear structure, V: Defense mechanisms, T: Signal transduction mechanisms, M: Cell wall/membrane biogenesis, N: Cell motility, Z: Cytoskeleton, W: Extracellular structures, U: Intracellular trafficking and secretion, O: Post-translational modification, protein turnover, chaperones, C: Energy production and conversion, G: Carbohydrate transport and metabolism, E: Amino acid transport and metabolism, F: Nucleotide transport and metabolism, H: Coenzyme transport and metabolism, I: Lipid transport and metabolism, P: Inorganic ion transport and metabolism, Q: Secondary metabolites biosynthesis, transport, and catabolism, R: General function prediction only, S: Function unknown.

(p-value = 0.001). A lower G + C content was correlated with the presence of an efflux pump.

The relationship between genome size and G + C content was examined in these strains, and no significant correlation was found.

## Quantitative principal component analysis

This preliminary exploratory analysis suggested non-random associations between genotypes and plasmid types, animal or human host, geography, and clinical symptoms (Figure 5). For example, MST 60, which corresponds to the plasmid QpDG, has only been isolated once from a rodent but never from humans. Plasmidless strains were associated with MST21 and Canada.

## Logistic (exponential) principal component analysis

We observed that the dataset mostly consisted of endocarditis cases from France (most strains were isolated from cardiac valves from French patients). The population of strains was organized into four groups identified by plasmids. QpDG was associated with animals such as rodents, cats, dogs, and MST60 (Figure 6). QpDV and MST1 were associated with acute Q fever (Figure 6). The association between QpDV and acute Q fever was also definitively confirmed in Tables 3, 4 and Figures 6–8. The QpH1 plasmid was rather associated with endocarditis in France. QpRS had an opposite behavior and was associated with the absence of acute Q fever, and this was

consistent among all analyses (Tables 3, 4 and Figures 5–8). Pericarditis, pneumonia, hepatitis, and fever were associated with acute Q fever (Figures 5, 6). Cat, dog, goat, mouse, and rodent strains grouped together (Figures 5, 6). Strains from Canada and the USA were mainly animal strains. These two PCAs (Pearson and logistic) prompted us to formally confirm an association between genotype and geography and clinical presentation. Statistical association with the host was not further tested because most strains [178/205 (86.8%)] were isolated from humans.

## Genotypes and geography

Among the strains with country of origin and MST information, we were able to identify many significant associations between genotypes and countries of origin (Supplementary Table 2). The strongest associations observed were between MST21 (without plasmid) and Canada (n = 5, Pearson coefficient 0.84, *p* < 0.05), between MST33 and The Netherlands (n = 10, 0.67, *p* < 0.05), and between MST17 and French Guyana (n = 3, 0.57, *p* < 0.05), thus confirming the existence of *C. burnetii* genotypes (Table 5).

## Genotypes and acute versus persistent focalized infections

Among 173 human strains for which the plasmid type and clinical presentation was available, we found a significant difference between strains. Strains with a QpRS plasmid [0/55 (0%)] and strains without plasmid [0/3 (0%)] were never
<figure><img src="image_5.png" /><figcaption>B</figcaption></figure>

<figure><img src="image_6.png" /><figcaption>C</figcaption></figure>

<figure><img src="image_7.png" /><figcaption>FIGURE 3
Gene variation. (A) Core/pangenome ratio (B): total genes; (C) core gene variation with the addition of new genome sequences.</figcaption></figure>

isolated from patients with acute Q fever. QpDV and QpH1 plasmid-positive strains were isolated from both patients with acute Q fever and persistent focalized infections. Acute Q fever corresponded to 17/44 (38.6%) QpDV strains but only 13/71 (18.3%) QpH1 strains. The difference in proportion was very significant between the four types (Chi-square = 26.1; *p* < 0.0001). When analyzing case significance with the Fisher's exact test, non-plasmid strains and QpRS strains were significantly associated with persistent focalized infections

(*p* < 0.05); QpRS strains were significantly associated with the absence of acute Q fever; QpDV strains were significantly associated with more acute Q fever and fewer persistent focalized infections (*p* < 0.05–Figures 7, 8).

After selecting the human strains for which the plasmid and clinical presentation were available, we excluded plasmids represented by less than four strains (plasmidless strains, *n* = 3) and clinical forms represented by less than four strains (pneumonia, lymphatic infection). We therefore
<figure><img src="image_8.png" /></figure>
<figure><img src="image_9.png" /><figcaption>FIGURE 5
Principal component analysis with genomic information [genome size, G + C content, multispacer sequence typing (MST), and plasmid type], geography, host, and clinical presentation.</figcaption></figure>

included 158 strains in this analysis. The distribution between plasmids (QpH1, QpDV, and QpRS) and clinical presentations (acute Q fever, endocarditis, vascular infection, miscarriage, osteoarticular infection) was described in **Figure 7** and **Table 6**.

Using the chi-square test, the Fisher's exact test per case (Table 3) and residual analysis [Table 4 and Figure 8; XLSTAT 2019 (Addinsoft, Paris, France)], we found that clinical presentation and plasmid type were not independent (Observed Chi-square 38.4, critical value 15.5, *p* < 0.0001). Strains with a QpDV plasmid were associated with more acute Q fever, less endocarditis and less osteoarticular infection. In addition, QpRS strains were not associated with any acute Q fever (0/53) but were associated with endocarditis. There was neither significant association for the other clinical forms nor for QpH1 with any clinical presentation. The QpRS plasmid was never identified in any acute Q fever strain but no specific geographic association was identified as the 63 QpRS plasmid-positive strains in our study were identified in several countries (49 in France, six in the USA, two in Spain, two in Russia, one in French Guiana, one in Senegal, one in the United Kingdom, one in Uzbekistan).

Discussion

We performed a pangenomic analysis of 75 *C. burnetii* isolates to understand better the evolutionary traits of their genomes. Focus was given on the links between genome content, epidemiological characteristics and pathogenicity. To our knowledge, the present work is the largest study on the association of *C. burnetii* strains and clinical presentation. Eighty-five percent of studied strains were isolated from human samples, including 69% from patients suffering from persistent focalized Q fever, and 59% from France. We first performed a phylogenomic analysis. Regardless of whether the studied genomes were complete or draft genomes, or whether they were sequenced using Illumina Miseq or Solid, clustering of the strains was consistent with the MST genotyping and the core genome SNP analysis, thus confirming the genotype-specific genetic differences found in other studies (Supplementary Figure 2; Kuley et al., 2017). The reliability of genotyping was also demonstrated during the Dutch outbreak (Kuley et al., 2017). Thus, the use of MST genotyping as a sensitive
<figure><img src="image_10.png" /><figcaption>FIGURE 6
Logistic (exponential) principal component analysis specifically adapted to the qualitative data.</figcaption></figure>

<table><caption>TABLE 3 Significance by cell.</caption><thead><tr><th>Plamsid</th><th>Acute Q fever</th><th>Endocarditis</th><th>Vascular infection</th><th>Osteoarticular infection</th><th>Miscarriage</th></tr></thead><tbody><tr><td>QpH1</td><td>></td><td>></td><td>&lt;</td><td>></td><td>&lt;</td></tr><tr><td>QpDV</td><td></td><td>&lt;</td><td>&lt;</td><td>&lt;</td><td>&gt;</td></tr><tr><td>QpRS</td><td>&lt;</td><td>&gt;</td><td>&gt;</td><td>&gt;</td><td>&lt;</td></tr></tbody></table>

Fisher's exact test. Values displayed in gray cells were significant at the level alpha = 0.05. XLSTAT2019 (Addinsoft, Paris, France).

<table><caption>TABLE 4 Residuals (adjusted).</caption><thead><tr><th></th><th>Acute Q fever</th><th>Endocarditis</th><th>Vascular infection</th><th>Osteoarticular infection</th><th>Miscarriage</th></tr></thead><tbody><tr><td>QpH1</td><td>0.036</td><td>0.65</td><td>-1.14</td><td>0.35</td><td>-0.32</td></tr><tr><td>QpDV</td><td>4.77</td><td>-4.06</td><td>-0.18</td><td>-1.38</td><td>1.82</td></tr><tr><td>QpRS</td><td>-4.32</td><td>2.95</td><td>1.36</td><td>0.87</td><td>-1.29</td></tr></tbody></table>

Values displayed in bold are significant at the level alpha = 0.05. XLSTAT2019 (Addinsoft, Paris, France).
<figure><img src="image_11.png" /><figcaption>FIGURE 7
Distribution between plasmids and clinical presentation.</figcaption></figure>

<figure><img src="image_12.png" /><figcaption>FIGURE 8
Plasmids and clinical residual analysis. Red: significant difference (p-value < 0.05).</figcaption></figure>
<table><caption>TABLE 5 Associations of multispacer sequence typing (MST) and countries.</caption><thead><tr><th>Country</th><th>Number of strains from this country</th><th>MST</th><th>Pearson's correlation coefficient</th></tr></thead><tbody><tr><td>Canada</td><td>5</td><td>21</td><td>0.84</td></tr><tr><td>The Netherlands</td><td>10</td><td>33</td><td>0.67</td></tr><tr><td>French Guiana</td><td>3</td><td>17</td><td>0.57</td></tr><tr><td>Russia</td><td>5</td><td>7</td><td>0.57</td></tr><tr><td>Russia</td><td>5</td><td>23</td><td>0.57</td></tr><tr><td>Senegal</td><td>4</td><td>19</td><td>0.50</td></tr><tr><td>Romania</td><td>3</td><td>18</td><td>0.42</td></tr><tr><td>USA</td><td>10</td><td>60</td><td>0.37</td></tr><tr><td>France</td><td>147</td><td>1</td><td>0.24</td></tr><tr><td>Senegal</td><td>4</td><td>6</td><td>0.20</td></tr><tr><td>USA</td><td>10</td><td>8</td><td>0.19</td></tr><tr><td>France</td><td>147</td><td>20</td><td>0.17</td></tr><tr><td>USA</td><td>10</td><td>16</td><td>0.17</td></tr><tr><td>Spain</td><td>3</td><td>4</td><td>0.16</td></tr><tr><td>Spain</td><td>3</td><td>8</td><td>0.16</td></tr></tbody></table>

Only significant association according to Pearson correlation test are shown (p < 0.05). Country with only one or two strains were excluded.

and less-time consuming tool should be considered in the future as a discriminatory tool for outbreak investigation and molecular epidemiology. In this study we identified a new MST genotype (MST60-Dugway strain) observed in a strain causing an independent cluster. In addition, we observed a concordance between MST genotyping, core genome analysis and plasmid type, which suggests an ancient acquisition of plasmids.

Comparative genomic analysis showed a high degree of genetic similarity between *C. burnetti* strains belonging to the same cluster. The highest values of OrthoANI and dDDH analyses were observed among strains from the same MST and harboring the same plasmid type. The classification of proteins by COG categories showed a homogenous distribution among many strains, especially for categories A, C, D, E, F, G, H, I, J, and K. Similarly, COG profiles of strains from a first phylogenetic cluster represented by four MST genotypes (MST12, MST13, MST32, MST33) exhibited greater gene numbers in the L-Z COG categories. These strains were isolated from Europe and had high G + C contents, numbers of genes and genome sizes than the remaining *Coxiella* genomes, thereby confirming close regional evolutionary relationships of genomic groups as shown

in previous studies, notably between UK strains and those from neighboring European countries (Hemsley et al., 2019).

Strain CB175, restricted to French Guiana, exhibited a COG profile as unique as its geographic distribution and virulence (Mahamat et al., 2013). Strain CB175 had the highest number of genes, irrespective of its smaller genome size compared to 67% of strains in this study. Thus, it was suspected that the virulence of this strain be linked to its genome reduction (D'Amato et al., 2015). The two strains CB94 and CB207, isolated from patients with hepatitis, exhibited similar COG profile (G-Z categories; Figure 2), the same MST16 genotype and the same G + C content despite a significant difference in genome size (Table 1). These observations indicate that the pathogenicity and clinical expression of *C. burnetti* strains may be related to their genomic content. The pangenome analysis also demonstrated a high genetic homogeneity among strains from the same MST type, which was in concordance with previous studies (D'Amato et al., 2014; Kuley et al., 2017). It also unveiled a low core/pangenome ratio, strain CB175 being an important contributor to this low value. This strain was demonstrated to have a 6,105-bp genomic deletion and a greater number of gene mutations when compared to strain NM-I (D'Amato et al., 2015). In addition, important core/pangenome percentage shifts were observed between MST groups and in some specific strains. By testing associations between strains characterized in the present study and including data from previous studies (Glazunova et al., 2005; Angelakis et al., 2013), we definitively demonstrated associations between genomic types (MST and plasmid types) and geographic distribution or clinical forms. Therefore MST genotyping may enable screening the strain type during an outbreak, which may play a key role in orientating the patient management and treatment approach. Although some genotypes appeared to be widely distributed, some were clearly specific, notably in French Guiana or during the Dutch outbreak. We also confirmed the association of plasmidless MST21 strains and Canada. This suggests that the geographic distribution may vary greatly among MST genotypes and that specific clonal complexes are associated with explosive outbreaks (MST33 in the Netherlands) or a localized hyperendemicity (MST17 in French Guiana) while other genotypes have broader geographical and temporal distributions, as is the case for MST1. Other studies also demonstrated a preferential associations of some genotypes with specific hosts (Roest et al., 2011) or specific epidemiological situations (Eldin et al., 2017).

<table><caption>TABLE 6 Association between plasmid type and clinical presentation.</caption><thead><tr><th>Plasmid</th><th>Acute Q fever</th><th>Endocarditis</th><th>Vascular infection</th><th>Osteoarticular infection</th><th>Miscarriage</th></tr></thead><tbody><tr><td>QpH1</td><td>13</td><td>45</td><td>4</td><td>3</td><td>3</td></tr><tr><td>QpDV</td><td>17</td><td>13</td><td>3</td><td>0</td><td>4</td></tr><tr><td>QpRS</td><td>0</td><td>42</td><td>7</td><td>3</td><td>1</td></tr></tbody></table>
Previous studies showed a relation between the presence of QpH1 and QpDV plasmids and acute infections or abortions (Angelakis et al., 2013). In our study, we observed that all isolates causing acute infections harbored either the QpH1 or QpDV plasmids and never the QpRs plasmid, even though these plasmid types are not exclusively found in acute Q fever-associated strains (Eldin et al., 2017). This finding has a clinical impact because our results suggest that QpRS-harboring *C. burnetii* strains cause asymptomatic primary infections. If this is confirmed, it would be impossible to prevent endocarditis or other persistent focused infections after primary infections in patients infected by such strains. This finding was not related to the confounding role of geography as several QpRS strains are found in several countries. Based on these findings, confirmed by genomic and statistical analyses, we assume that strain characteristics play an important role in the clinical manifestations of Q fever.

As previously described, *C. burnetii* possesses an important number of IS sequences in comparison to other intracellular bacteria (Seshadri et al., 2003). We detected IS sequences in all studied genomes which is in concordance with other studies that observed numerous IS sequences in *C. burnetii* genomes, thus conferring a greater genome plasticity to *C. burnetii* in comparison to other intracellular bacteria (Beare et al., 2009; Kuley et al., 2017). Genomic plasticity was clearly demonstrated in Figure 4. The detection of genes present in other intracellular bacteria, as is the case of the unique genes found in strains CB5, CB193, and Dugway, suggests that, during its evolution, *C. burnetii* shared the same host with other intracellular bacteria, thus enabling gene transfers.

Regarding antibiotic susceptibility, rare resistance issues are observed in clinical practice. Genes coding for resistance to betalactams and fluoroquinolones were found in all strains (Spyridaki, 2002; Vranakis et al., 2010). The presence of an efflux pump was observed in some strains that belong to the same branch in the phylogenetic tree. This efflux pump belonging to the MATE family was found in different bacteria such as *Escherichia coli* and *Klebsiella pneumoniae* but has not been associated with any antibiotic resistance to date (Ogawa et al., 2015). Statistical analyses showed a significant correlation between the presence of this efflux pump and a low G + C content, suggesting that it might be acquired via horizontal gene transfer, thus reinforcing the idea of genome plasticity and genotype-specific pathogenicity.

In conclusion, genomic analysis enabled us to understand more of the genomic dynamics of *C. burnetii*, its evolution, host adaptation, epidemiology, and clinical profile. In this study we demonstrated that *C. burnetii* exhibits a high genome plasticity and that some of the genomic features are correlated to the clinical outcome, epidemiological characteristics and geographic distribution. Future studies including more animal and human-infecting strains from different MST genotypes and from additional geographical locations will complete

our current perspective of *C. burnetii* genomic evolution and pathogenesis.

## Data availability statement

All sequences presented in this study were deposited in GenBank under the accession numbers listed in Table 1.

## Author contributions

RA and P-EF designed the project. RA, HA, JD, AD, AC, ER, KS, and RZ performed the experiments and genome analysis. RA and MM performed the statistical analysis. RA performed the data analysis. RA, MM, and P-EF wrote the manuscript. P-EF, DR, and MM revised the manuscript. All authors contributed to the article and approved the submitted version.

## Funding

This study was supported by the Méditerranée Infection foundation, the National Research Agency under the program "Investissements d'avenir," reference ANR-10-IAHU-03 and by Région Provence Alpes Côte d'Azur and European funding FEDER PRIMI.

## Conflict of interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

## Publisher's note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

## Supplementary material

The Supplementary Material for this article can be found online at: https://www.frontiersin.org/articles/10.3389/fmicb.2022.1022356/full#supplementary-material
References

Angelakis, E., and Raoult, D. (2010). Q fever. *Vet. Microbiol.* 140, 297–309.
doi: 10.1016/j.vetmic.2009.07.016

Angelakis, E., Edouard, S., Lafranchi, M.-A., Pham, T., Lafforgue, P., and
Raoult, D. (2014). Emergence of Q fever arthritis in France. *J. Clin. Microbiol.* 52,
1064–1067. doi: 10.1128/JCM.03371-13

Angelakis, E., Million, M., D'Amato, F., Rouli, L., Richet, H., Stein, A.,
et al. (2013). Q fever and pregnancy: Disease, prevention, and strain specificity.
*Eur. J. Clin. Microbiol. Infect. Dis.* 32, 361–368. doi: 10.1007/s10096-012-1
750-3

Angelakis, E., Richet, H., Rolain, J. M., La Scola, B., and Raoult, D. (2012).
Comparison of real-time quantitative PCR and culture for the diagnosis of
emerging rickettsioses. *PLoS Negl. Trop. Dis.* 6:e1540. doi: 10.1371/journal.pntd.
0001540

Auch, A. F., von Jan, M., Klenk, H.-P., and Göker, M. (2010). Digital DNA-DNA hybridization for microbial species delineation by means of genome-to-genome sequence comparison. *Stand. Genom. Sci.* 2, 117–134. doi: 10.4056/sigs.53
1120

Bankevich, A., Nurk, S., Antipov, D., Gurevich, A. A., Dvorkin, M., Kulikov, A. S., et al. (2012). SPAdes - A new genome assembly algorithm and its applications to single-cell sequencing. *J. Comput. Biol.* 19, 455–477. doi: 10.1089/cmb.2012.0021

Beare, P. A., Unsworth, N., Andoh, M., Voth, D. E., Omsland, A., Gilk, S. D., et al. (2009). Comparative genomics reveal extensive transposon-mediated genomic plasticity and diversity among potential effector proteins within the genus *Coxiella*. *Infect. Immun.* 77, 642–656. doi: 10.1128/IAI.011
41-08

Bertelli, C., Laird, M. R., Williams, K. P., Simon Fraser University Research
Computing Group, Lau, B. Y., Hoad, G., et al. (2017). IslandViewer 4: Expanded
prediction of genomic islands for larger-scale datasets. *Nucleic Acids Res.* 45:W30-
W35. doi: 10.1093/nar/gkx343

Caron, F., Meurice, J. C., Ingrand, P., Bourgoin, A., Masson, P., Roblot, P., et al. (1998). Acute Q fever pneumonia: A review of 80 hospitalized patients. *Chest* 114, 808–813. doi: 10.1378/chest.114.3.808

Chen, L. (2004). VFDB: A reference database for bacterial virulence factors.
*Nucleic Acids Res.* 33:D325–D328. doi: 10.1093/nar/gki008

D'Amato, F., Eldin, C., and Raoult, D. (2016). The contribution of genomics to the study of Q fever. *Future Microbiol.* 11, 253–272. doi: 10.2217/fmb.1
5.137

D'Amato, F., Eldin, C., Georgiades, K., Edouard, S., Delerce, J., Labas, N., et al. (2015). Loss of TSS1 in hypervirulent *Coxiella burnetii* 175, the causative agent of Q fever in French Guiana. *Comp. Immunol. Microbiol. Infect. Dis.* 41, 35–41. doi: 10.1016/j.cimid.2015.04.003

D'Amato, F., Rouli, L., Edouard, S., Tyczka, J., Million, M., Robert, C., et al. (2014). The genome of *Coxiella burnetii* Z3055, a clone linked to the Netherlands Q fever outbreaks, provides evidence for the role of drift in the emergence of epidemic clones. *Comp. Immunol. Microbiol. Infect. Dis.* 37, 281–288. doi: 10.1016/j.cimid.2014.08.003

de Bruin, A., van Alphen, P. T. W., van der Plaats, R. Q. J., de Heer, L. N. D., Reusken, C. B. E. M., van Rotterdam, B. J., et al. (2012). Molecular typing of *Coxiella burnetii* from animal and environmental matrices during Q fever epidemics in the Netherlands. *BMC Vet. Res.* 8:165. doi: 10.1186/1746-6148-8-165

Derrick, E. H. (1983). "Q" fever, a new fever entity: Clinical features, diagnosis and laboratory investigation. *Clin. Infect. Dis.* 5, 790–800. doi: 10.1093/clinids/5.4.790

Eldin, C., Melenotte, C., Mediannikov, O., Ghigo, E., Million, M., Edouard, S., et al. (2017). From Q fever to *Coxiella burnetii* infection: A paradigm change. *Clin. Microbiol. Rev.* 30, 115–190. doi: 10.1128/CMR.00045-16

Fournier, P. E., Etienne, J., Harle, J. R., Habib, G., and Raoult, D. (2001). Myocarditis, a rare but severe manifestation of Q fever: Report of 8 cases and review of the literature. *Clin. Infect. Dis.* 32, 1440–1447. doi: 10.1086/320159

Frankel, D., Richet, H., Renvoisé, A., and Raoult, D. (2011). Q fever in France,
1985-2009. *Emerg. Infect. Dis.* 17, 350–356. doi: 10.3201/eid1703.100882

Glazunova, O., Roux, V., Freylikman, O., Sekeyova, Z., Fournous, G., Tyczka, J., et al. (2005). *Coxiella burnetii* genotyping. *Emerg. Infect. Dis.* 11, 1211–1217. doi: 10.3201/eid1108.041354

Hemsley, C. M., O'Neill, P. A., Essex-Lopresti, A., Norville, I. H., Atkins, T. P., and Titball, R. W. (2019). Extensive genome analysis of *Coxiella burnetii* reveals limited evolution within genomic groups. *BMC Genom.* 20:441. doi: 10.1186/s12864-019-5833-8

Hilbink, F., Penrose, M., Kovacova, E., and Kazar, J. (1993). Q-fever is absent from New-Zealand. *Int. J. Epidemiol.* 22, 945–949. doi: 10.1093/ije/22.5.945

Kuley, R., Bossers-deVries, R., Smith, H. E., Smits, M. A., Roest, H. I. J., and Bossers, A. (2015). Major differential gene regulation in *Coxiella burnetii* between *in vivo* and *in vitro* cultivation models. *BMC Genom.* 16:953. doi: 10.1186/s12864-015-2143-7

Kuley, R., Kuijt, E., Smits, M. A., Roest, H. I. J., Smith, H. E., and Bossers, A. (2017). Genome plasticity and polymorphisms in critical genes correlate with increased virulence of dutch outbreak-related *Coxiella burnetii* strains. *Front. Microbiol.* 8:1155. doi: 10.3389/fmicb.2017.01526

Lechner, M., Findeiß, S., Steiner, L., Marz, M., Stadler, P. F., and Prohaska, S. J. (2011). Proteinortho : Detection of (Co-)orthologs in large-scale analysis. *BMC Bioinform.* 12:124. doi: 10.1186/1471-2105-12-124

Lee, I., Kim, Y. O., Park, S.-C., and Chun, J. (2016). OrthoANI: An improved algorithm and software for calculating average nucleotide identity. *Int. J. Syst. Evol. Microbiol.* 66, 1100–1103. doi: 10.1099/ijsem.0.000760

Lyytikäinen, O., Ziese, T., Schwartländer, B., Matzdorff, P., Kuhnen, C., Jäger, C., et al. (1998). An outbreak of sheep-associated Q fever in a rural community in Germany. *Eur. J. Epidemiol.* 14, 193–199. doi: 10.1023/A:1007452503863

Mahamat, A., Edouard, S., Demar, M., Abboud, P., Patrice, J.-Y., La Scola, B., et al. (2013). Unique clone of *Coxiella burnetii* causing severe Q fever, French Guiana. *Emerg. Infect. Dis.* 19, 1102–1104. doi: 10.3201/eid1907.130044

Manfredi Selvaggi, T., Rezza, G., Scagnelli, M., Rigoli, R., Rassu, M., De Lalla, F., et al. (1996). Investigation of a Q-fever outbreak in northern Italy. *Eur. J. Epidemiol.* 12, 403–408. doi: 10.1007/BF00145305

Melenotte, C., Million, M., and Raoult, D. (2020). New insights in *Coxiella burnetii* infection: Diagnosis and therapeutic update. *Expert Rev. Anti Infect. Ther.* 18, 75–86. doi: 10.1080/14787210.2020.1699055

Melenotte, C., Million, M., Audoly, G., Gorse, A., Dutronc, H., Roland, G., et al. (2016). B-cell non-Hodgkin lymphoma linked to *Coxiella burnetii*. *Blood* 127, 113–121. doi: 10.1182/blood-2015-04-639617

Million, M., and Raoult, D. (2017). No such thing as chronic Q fever. *Emerg. Infect. Dis.* 23, 856–857. doi: 10.3201/eid2305.151159

O'Donnell, M. E., Manshani, N., McCaughey, C., Soong, C., and Lee, B. (2007).
*Coxiella burnetii* infection of an aortic graft with multiple vertebral body erosion.
*J. Vasc. Surg.* 45, 399–403. doi: 10.1016/j.jvs.2006.09.016

Ogawa, W., Minato, Y., Dodan, H., Onishi, M., Tsuchiya, T., and Kuroda, T. (2015). Characterization of Mate-type multidrug efflux pumps from *Klebsiella pneumoniae* MGH78578. *PLoS One* 10:e0121619. doi: 10.1371/journal.pone.0121619

Page, A. J., Taylor, B., Delaney, A. J., Soares, J., Seemann, T., Keane, J. A., et al. (2016). SNP-sites: Rapid efficient extraction of SNPs from multi-FASTA alignments. *Microb. Genom.* 2:e000056. doi: 10.1099/mgen.0.000056

Raoult, D. (2017). Q Fever: Confusion between chronic infection and chronic fatigue. *Clin. Infect. Dis.* 65, 1054–1055. doi: 10.1093/cid/cix469

Raoult, D., Marrie, T. J., and Mege, J. L. (2005). Natural history and pathophysiology of Q fever. *Lancet Infect. Dis.* 5, 219–226. doi: 10.1016/S1473-3099(05)70052-9

Roest, H. I. J., Ruuls, R. C., Tilburg, J. J. H. C., Nabuurs-Franssen, M. H., Klaassen, C. H. W., Vellema, P., et al. (2011). Molecular epidemiology of *Coxiella burnetii* from ruminants in Q fever outbreak, the Netherlands. *Emerg. Infect. Dis.* 17, 668–675. doi: 10.3201/eid1704.101562

Samuel, J. E., Frazier, M. E., and Mallavia, L. P. (1985). Correlation of plasmid type and disease caused by *Coxiella burnetii*. *Infect. Immun.* 49, 775–779. doi: 10.1128/iai.49.3.775-779.1985

Seemann, T. (2014). Prokka: Rapid prokaryotic genome annotation.
*Bioinformatics* 30, 2068–2069. doi: 10.1093/bioinformatics/btu153

Seshadri, R., Paulsen, I. T., Eisen, J. A., Read, T. D., Nelson, K. E., Nelson, W. C., et al. (2003). Complete genome sequence of the Q-fever pathogen
Coxiella burnetii. Proc. Natl. Acad. Sci. U.S.A. 100, 5455–5460. doi: 10.1073/pnas.0931379100

Sidi-Boumedine, K., Ellis, R. J., Adam, G., Prigent, M., Angen, O., Aspán, A., et al. (2014). Draft genome sequences of six ruminant *Coxiella burnetii* isolates of European origin. *Genome Announc.* 2:e00285-14. doi: 10.1128/genomeA.00285-14

Smith, D. L., Ayres, J. G., Blair, I., Burge, P. S., Carpenter, M. J., Caul, E. O., et al. (1993). A large Q fever outbreak in the West Midlands: Clinical aspects. *Respir. Med.* 87, 509–516. doi: 10.1016/0954-6111(93)90006-L

Spyridaki, I. (2002). Mechanisms of resistance to fluoroquinolones in *Coxiella burnetii*. *J. Antimicrob. Chemother.* 49, 379–382. doi: 10.1093/jac/49.2.379

Tettelin, H., Riley, D., Cattuto, C., and Medini, D. (2008). Comparative genomics: The bacterial pan-genome. *Curr. Opin. Microbiol.* 11, 472–477. doi: 10.1016/j.mib.2008.09.006

Vranakis, I., Sandalakis, V., Chochlakis, D., Tselentis, Y., and Psaroulaki, A. (2010). DNA gyrase and topoisomerase IV mutations in an *in vitro* fluoroquinolone-resistant *Coxiella burnetii* strain. *Microb. Drug Resist.* 16, 111–117. doi: 10.1089/mdr.2010.0015