Legionella polyplacis genomes reconstruction and comparison
The contigs corresponding to Legionella polyplacis were visually identified using ORF prediction done in the Geneious package (the prokaryotic gene arrangement could be readily recognized by the density of predicted ORFs). In each sample, except for the DPH41, the genome ofL. polyplacis was assembled into a single complete contig. In the DPH41 assembly, we were not able to find the symbiont genome, despite the good assembly quality of the louse genome. Complete L. polyplacis genomes were annotated in RAST (Aziz et al., 2008) and aligned using Mafft algorithm implemented in Geneious. In three samples (Ne125b_Kot_SW, 99b_Pro_SE, 98c_Pro_SE) the rRNA region did not assemble correctly and was only represented by short fragments. To extend these fragments into full length, we used the program aTRAM 2.0 (Allen, LaFrance, Folk, Johnson, & Guralnick, 2018). Phylogenetic analysis of the resulting 530,042 bp long matrix was performed in Phyml (Guindon et al., 2010). The evolutionary model GTR+I+G was determined for the whole concatenated matrix (considering the strong phylogenetic signal in the matrix, the model selection is a purely formalistic step in the presented analysis) by a corrected Akaike information criterion in jModelTest2 (Darriba et al., 2012) based on the AIC. Comparison of gene content was done manually using the annotated alignment. We were taking into consideration the following criteria: presence of the annotations across all genomes; presence of the corresponding sequence (in case of missing annotation); distribution of differences (i.e. does the difference reflect the SE/SW split).