Host adaptation and microbial competition drive Ralstonia solanacearum phylotype I evolution in the Republic of Korea.

Bacterial wilt caused by the Ralstonia solanacearum species complex (RSSC) threatens the cultivation of important crops worldwide. We sequenced 30 RSSC phylotype I (R. pseudosolanacearum) strains isolated from pepper (Capsicum annuum) and tomato (Solanum lycopersicum) across the Republic of Korea. These isolates span the diversity of phylotype I, have extensive effector repertoires and are subject to frequent recombination. Recombination hotspots among South Korean phylotype I isolates include multiple predicted contact-dependent inhibition loci, suggesting that microbial competition plays a significant role in Ralstonia evolution. Rapid diversification of secreted effectors presents challenges for the development of disease-resistant plant varieties. We identified potential targets for disease resistance breeding by testing for allele-specific host recognition of T3Es present among South Korean phyloype I isolates. The integration of pathogen population genomics and molecular plant pathology contributes to the development of location-specific disease control and development of plant cultivars with durable resistance to relevant threats.

Data Summary

The genome assemblies generated in this study have been deposited in the National Center for Biotechnology Information (NCBI) database under BioProject number PRJNA593908; accession numbers are listed in Table 1. The other publicly available assemblies included in this work are listed in Table S1 (available in the online version of this article).

Table 1.

phylotype I isolates and assembly statistics

Strain	Year	Location of isolation	Host of isolation	Assembly length (bp)	Scaffolds†	Longest scaffold	N50	Coverage	Accession
Pe_1*	2000	Chungnam Seosan, Republic of Korea	C. annuum	5 839 983	5	3 681 740	3 681 758	306	WSNR
Pe_2	2000	Chungbuk Cheongwon, Republic of Korea	C. annuum	5 600 524	213	254 641	85 236	289	WSYW
Pe_3*	2001	Chungbuk Chungju, Republic of Korea	C. annuum	5 938 020	3	3 851 714	3 851 698	500	WSYS
Pe_4	2001	Jeju Bukjeju, Republic of Korea	C. annuum	5 593 270	203	317 826	77 130	247	WSYZ
Pe_9	2002	Chungnam Cheongyang, Republic of Korea	C. annuum	5 723 071	278	233 723	78 017	240	WSYY
Pe_13	2002	Jeonbuk Imsil, Republic of Korea	C. annuum	5 865 496	248	214 411	77 725	208	WSYX
Pe_15	2002	Jeonbuk Imsil, Republic of Korea	C. annuum	6 033 637	237	329 719	73 852	291	WSZA
Pe_18	2002	Jeonnam Haenam, Republic of Korea	C. annuum	6 029 178	254	241 872	83 596	261	JAABKA
Pe_24	2002	Chungnam Gongju, Republic of Korea	C. annuum	5 917 705	253	421 926	85 211	319	JAABJZ
Pe_26	2003	Chungnam Seosan, Republic of Korea	C. annuum	5 747 672	200	235 891	77 520	286	JAABJY
Pe_27*	2003	Chungbuk Cheongwon, Republic of Korea	C. annuum	5 853 871	9	1 965 711	1 655 033	450	WSYT
Pe_28	1999	Gyeonggi Hwaseong, Republic of Korea	C. annuum	5 830 258	247	236 369	74 138	243	JAABJX
Pe_30	2005	Jeonbuk Imsil, Republic of Korea	C. annuum	5 894 388	227	246 051	83 264	271	JAABJW
Pe_39*	2005	Jeonbuk Jeongeup, Republic of Korea	C. annuum	5 919856	4	3 808 376	3 808 376	490	WSYU
Pe_42	2005	Jeonnam Naju, Republic of Korea	C. annuum	5 788 311	272	170 277	74 548	281	JAABJV
Pe_45	2005	Chungnam Gongju, Republic of Korea	C. annuum	5 702 018	211	317 902	72 734	246	JAABJU
Pe_49	2005	Chungnam Cheongyang, Republic of Korea	C. annuum	5 645 699	223	241 872	73 220	300	JAABJT
Pe_51	2005	Chungnam Taean, Republic of Korea	C. annuum	5 722 228	229	406 488	86 970	285	JAABJS
Pe_52	2005	Chungbuk Goesan, Republic of Korea	C. annuum	5 688 143	240	205 520	78 443	250	JAABJR
Pe_56	2005	Chungbuk Eumseong, Republic of Korea	C. annuum	5 642 649	225	292 052	86 052	258	JAABJQ
Pe_57*	2005	Gyeongbuk Cheongsong, Republic of Korea	C. annuum	5 720 795	2	3 754 354	3 754 354	490	WSYV
Pe_61	2005	Chungbuk Goesan, Republic of Korea	C. annuum	5 719 158	246	269 792	82 771	247	JAABJP
To_1	2008	Gangwon Hoengseong, Republic of Korea	S. lycopersicum	5 878 275	235	216 471	84 441	245	JAABJO
To_7	2008	Gangwon Hwacheon, Republic of Korea	S. lycopersicum	5 666 531	228	269 904	78 595	215	JAABJN
To_22	2008	Gangwon Cheorwon, Republic of Korea	S. lycopersicum	5 357 968	223	323 913	62 467	232	JAABJM
To_28	2008	Gangwon Yanggu, Republic of Korea	S. lycopersicum	5 656 151	235	275 581	69 771	245	JAABJL
To_36	2008	Gangwon Pyeongchang, Republic of Korea	S. lycopersicum	5 785 858	252	213 926	79 781	238	JAABJK
To_42	2008	Gangwon Chuncheon, Republic of Korea	S. lycopersicum	5 811 979	306	214411	61 263	112	JAABJJ
To_53	2008	Gangwon Hongcheon, Republic of Korea	S. lycopersicum	5 357 324	231	297 172	54 396	226	JAABJI
To_63	2008	Gyeongbuk Bonghwa, Republic of Korea	S. lycopersicum	5 612 635	208	210 182	82 829	203	JAABJH

*Hybrid assemblies.

†Contigs shorter than 500 bp not included.

The species complex causes frequent outbreaks of bacterial wilt disease in staple crops around the world. We found that diverse lineages of phylotype I are circulating across the Republic of Korea. These isolates have extensive repertoires of type III secreted effectors (T3Es) and are subject to frequent recombination. Recent recombination hotspots overlap with loci encoding contact-dependent inhibition systems, which mediate interbacterial antagonism. Our examination of T3E diversity across the entire species complex indicates that some effectors are subject to frequent horizontal transfer. Selection imposed by host recognition may result in the loss or diversification of recognized pathogen effectors. Indeed, we found some T3Es among the South Korean isolates elicit allele-specific recognition responses in planta. This work provides an insight into how host adaptation and microbial competition may both drive evolution.

Introduction

is a soil-borne pathogen that causes bacterial wilt disease of over 200 species in 50 different plant families, presenting a major threat to global crop production. This β-proteobacterium colonizes vascular tissues by invading root systems via wounds. Rapid proliferation and biofilm formation block water flow, causing wilt and death [1, 2]. causes severe outbreaks in crops of critical economic and nutritional importance, including banana (Musa spp.), potato (Solanum tuberosum), tomato (S. lycopersicum) and pepper (Capsicum annuum) [1]. The broad host range, extensive distribution and genetic heterogeneity of present major obstacles for the design of sustainable disease management strategies.

was initially described as a single species within prior to reassignment to [3, 4]. Low (<70 %) levels of DNA hybridization between strains indicated that is more appropriately described as a species complex (RSSC) [5]. Four RSSC phylotypes were delineated based on 16S–23S rRNA gene intergenic spacer (ITS) regions and the hrpB and egl genes, coupled with comparative genomic hybridization analysis. Phylotype I is mainly present in Asia, phylotype II mainly in the Americas, phylotype III in Africa and phylotype IV in Japan, Indonesia, Australia and the Philippines [1]. Multi-locus sequence analysis and whole-genome comparisons indicate that phylotype I and III are more closely related to each other, while phylotype IV is most divergent. Although we refer to phylotype designations in this paper, taxonomic revisions have been proposed to reclassify the RSSC into three separate species: (phylotypes I and III), (phylotype II) and (phylotype IV) [6, 7].

phylotypes have ~5 Mb bipartite genomes comprising two replicons 3.5–3.7 and ~1.6–2.1 Mb in size. The smaller replicon is conventionally referred to as a megaplasmid, although both replicons encode essential and pathogenesis-related genes. Multipartite genomes are frequently observed among genera harbouring soil and marine bacteria interacting with eukaryotic hosts [5]. The mosaic structure of RSSC genomes attests to the importance of recombination in enhancing the genetic variation evident in this species complex [5, 8]. Natural competence is considered a ubiquitous trait among RSSC phylotypes and has been demonstrated both in vitro and in planta [9, 10]. Coupat et al. [10] found that 80 % of 55 RSSC strains tested developed competence in vitro and exhibited single or multiple homologous recombination events up to 90 kb in length. High rates of recombination enhance the potential for rapid adaptation and the emergence of novel pathogen variants [8, 11].

Virulence mechanisms in the RSSC include biofilm production, secretion of cell wall-degrading enzymes, and type III secreted virulence effector (T3E) proteins. Coordinated synthesis of T3Es is achieved by the activity of HrpG/HrpB transcriptional regulators with respect to the hypersensitive response and pathogenicity (hrp) box cis-element in the promoter of T3E genes. RSSC T3Es target various compartments and functions to suppress host resistance and promote bacterial proliferation [12]. The RSSC has a vastly expanded type III effector repertoire compared to other well-studied plant pathogenic bacteria, with over 100 T3E families identified thus far [13]. In contrast, 70 and 53 orthologous effector families are reported in and spp., respectively [14, 15].

Plants have evolved intracellular receptors that monitor for the presence of pathogen effectors and activate a robust defence response, often culminating in a form of programmed cell death (referred to as a hypersensitive response) around the site of infection and recognition [16]. Effector-dependent host-specific interactions have been identified in R. solanacearum. RipAA (formerly avrA) was the first avirulence locus identified in NC252 triggering recognition in tobacco plants [17]. Further studies revealed that RipAA, RipP1 (formerly popP1) and RipB in the phylotype I strain GMI1000 are recognized in Nicotiana tabacum and Nicotiana benthamiana [18, 19]. A comprehensive survey of phylotype II BS048 effector recognition identified RipAA, RipE1 and RipH2 as cell death inducers in Nicotiana spp., tomato and lettuce, respectively [20]. The phylotype I GMI1000 RipE1 allele also triggers hypersensitive response in N. benthamiana [21, 22]. RipA family AWR motif-containing effectors are strong inducers of cell death in Nicotiana spp. [23]. The loss of RipAA was recently found to be correlated with the emergence of phylotype IIB strains pathogenic on cucurbits; related strains retaining RipAA are pathogenic on banana but not cucurbits [24].

We sequenced 30 RSSC phylotype I isolates from tomato and capsicum to identify the population structure and virulence factor repertoire of pathogens limiting production of these key crops in the Republic of Korea. Location and host of isolation do not have a strong impact on the population structure of South Korean phylotype I isolates; diverse phylotype I isolates capable of infecting solanaceous crops are broadly distributed across the Republic of Korea. Loci linked with competition are frequent targets of recombination, suggesting that microbial interactions play an important role in evolution. We examined South Korean phylotype I virulence factor evolution in the context of the global pool of RSSC T3Es and found that most strains carry large effector repertoires with evidence of frequent horizontal transfer. We identified T3Es that elicit stable host recognition and defence responses, regardless of which allelic variant is used, indicating they would be suitable targets for disease resistance breeding efforts. We also identified allelic variation in two T3Es resulting in the loss of host recognition. Rapid diversification of secreted effectors present challenges for the development of resistant plant varieties. The integration of pathogen population genomics with molecular plant pathology contributes to location-specific pathogen control strategies, as well as the development of cultivars with durable resistance to relevant threats.

Methods

Genome assembly and annotation

(Rso) were isolated from field-grown pepper and greenhouse tomato plants between 1999 and 2008 (Table 1). Tomato isolates were collected from the north in 2008, while pepper isolates were collected from central and southern regions of the South Korean peninsula. One isolate (Pe_4) was collected from Jeju island, southeast of the peninsula. Genomic DNA was extracted using the Promega Wizard Genomic DNA Purification kit and 150 bp paired-end sequencing was run on the Illumina NextSeq platform using Nextera and the library preparation protocol described by Baym et. al. [25]. Reads were filtered to remove adapter sequences (parameters: ktrim=r, k=23, mink=11, hdist=1, minlen=50, ftm=5 tpe tbo) and common sequence contaminants according to the National Center for Biotechnology Information (NCBI) Univec database (parameters: k=31, hdist=1), and quality-trimmed using the bbduk program available in the BBmap package (parameters: qtrim=r, trimq=10, minlen=50, maq=10 tpe).

Filtered and quality-trimmed reads were used for genome assembly with SPAdes v3.10.0, using --careful assembly with kmer size distribution –k 21, 33, 55 and 77 [26]. The coverage cutoff was set to auto. Coverage was assessed by mapping reads to the assembly using Bowtie 2, removing duplicates [27]. Assembly improvement was performed using Pilon v1.14 and assembly statistics were obtained using assembly-stats (https://github.com/sanger-pathogens/assembly-stats) [28]. Genomes were annotated with Prokka v.1.12 beta using a custom annotation database created using all 71 annotated assemblies available on GenBank 11 November 2017 [29]. Sixty-eight publicly available genomes were also reannotated with Prokka v1.13.3 for inclusion in downstream analyses (Table S1).

Nanopore sequencing was performed for five isolates: Pe_1, Pe_3, Pe_27, Pe_39 and Pe_57. DNA extraction was performed using Promega Wizard, and barcoding and library preparation were performed using the Native Barcoding Expansion 1–12 (EXP-NBD104) kit, in conjugation with the Ligation Sequencing kit (SQK-LSK109), according to the manufacturer’s recommendations. After a sequencing run of 48 h, basecalling was performed using Guppy followed by hybrid assembly with Unicycler using Illumina reads for each respective strain [30]. Genomes were annotated with Prokka v1.13.3 using the custom annotation database.

Core genome phylogeny of Korean phylotype I isolates

The core genome phylogeny of South Korean phylotype I isolates was obtained by mapping Illumina reads to the completed nanopore assembly of Pe_57. Regions of the reference genome corresponding to mobile genetic elements were removed after identification using the PHASTER phage and ISFinder databases; no integrative and conjugative elements were identified with ICEberg [31-33].

Readmapping of filtered, quality-trimmed reads was performed using BWA mem 0.7.17 r1188 [34]. Samtools v1.3.1 was used to sort and index reads, remove duplicate reads and generate pileup data [31]. Variants were called with BCFtools v1.3.1 (bcftools mpileup --consensus-caller --ploidy 1 --min-MQ 60) and a set of filtered single-nucleotide polymorphisms (SNPs) obtained using BCFTools v0.1.13 to retain SNPs with a minimum read depth of 10 in both directions, non-reference allele frequency of 0.95, minimum quality of 30 and minimum distance from a gap of 10 bp [35]. BEDtools v2.26.0 was used to determine positions where the depth of coverage fell below 10 reads [36]. A final consensus relative to the reference was obtained with BCFtools consensus v1.3.1, masking low-coverage positions identified by BEDtools, retaining both invariant sites and SNPs passing coverage, non-reference allele frequency and quality thresholds. An alignment was generated using the consensus sequences of each isolate. Alignment positions with one or more gaps or overlapping mobile genetic elements in the reference genome were removed. Initial tree building was performed with RAxML v7.2.8 using the parameters: -f a -p $RANDOM -x $RANDOM -# 100 m GTRGAMMA -T 8 [37].

ClonalFrameML v1.178 was employed to identify polymorphisms introduced via homologous recombination using the gap-free core genome alignment and initial maximum-likelihood phylogeny as input [38]. Alignment regions identified as likely to have been introduced via recombination were removed from the alignment and tree building was repeated using the non-recombinant alignment in RAxML v7.2.8 (parameters as above). Trees were visualized using FigTree v1.4.3, displaying only nodes with bootstrap support scores of 70 and above (https://github.com/rambaut/figtree/). This approach was employed to obtain a core genome phylogeny and alignment using the primary and secondary replicons of R. solanacaerum Pe_57.

Recent and ancestral recombination events among South Korean strains phylotype I strains were visualized using maskrc-svg.py with the core genome alignment along with the recombinant regions predicted by ClonalFrameML (https://github.com/kwongj/maskrc-svg). Recombination events are referred to as ancestral when predicted to occur at internal (ancestral) nodes, and extant or recent events when impacting on terminal nodes (individual strains). The number of recent recombinant events was calculated with pyGenomeTracks by counting the number of extant (recent) recombination events overlapping with annotated coding sequences in the reference genome GMI1000 chromosome and megaplasmid, excluding recent events in GMI1000 and Pe_57 (https://github.com/deeptools/pyGenomeTracks).

species complex phylogeny

The species complex tree was generated using an alignment of concatenated single-copy core genes. All Korean isolates and 68 publicly available genomes were annotated with Prokka (as above). OrthoFinder was implemented to identify single-copy orthologues [39]. Single-copy orthologue nucleotide sequences were employed for codon-aware alignment with PRANK [40]. Core gene alignments were concatenated and stripped of gap and invariant positions. Maximum-likelihood phylogeny reconstruction was then performed with RAxML v8.2.10 (parameters: -f a -p $RANDOM -x $RANDOM -# 100 m GTRCATX -T 16) [36]. Trees were visualized using FigTree v1.4.3, only displaying nodes with bootstrap support scores of 70 and above (https://github.com/rambaut/figtree/). Estimation of the population structure was performed with RhierBAPS using the core SNP alignment with two levels of clustering and setting the number of initial clusters to 20 [41].

Type III effector identification and comparative analyses

Effector predictions were performed on the South Korean phylotype I isolates and 68 publicly available genomes. First, sequences from the manually curated effector (RalstoT3E) database were downloaded from http://iant.toulouse.inra.fr/bacteria/annotation/site/prj/T3Ev3 and used as queries in tblastn searches of the 30 Korean genomes, retaining all hits with a minimum coverage of 30 % and pairwise identity of 50 %. Overlapping hits on the subject sequence were merged to generate the longest hit using BEDtools. The merged hit was then scanned to identify the predicted open reading frame (ORF), searching for start and stop codons and the presence of any premature stop codons. The presence of putative hrp boxes up to 500 bp upstream of the effector candidate ORF start codon was retained. Merged hits containing no start codons and/or premature stop codons were adjusted to reconstitute a functional ORF if possible. Start codon adjustment was performed by searching between the merged hit and putative hrp box (or up to 500 bp upstream of the merged hit for a start codon). Effector candidates were classified as pseudogenes in the absence of a start codon or in the presence of a stop codon resulting in a truncated ORF 30 % shorter than the shortest query sequence generating the initial hit.

Individual effector phylogenies were generated using PRANK codon-aware alignments followed by tree building using RAxML 8.2.12. Gene tree topologies were compared with the species complex core gene phylogeny using the normalized Robinson–Foulds metric implemented in ETE3, collapsing identical sequences into the same leaf and disregarding comparison between nodes with bootstrap scores lower than 70 [42].

Bacterial strains and plant growth conditions

DH5α and AGL1 were grown in Luria–Bertani (LB) medium supplied with appropriate antibiotics at 37 and 28 °C, respectively. N. benthamiana and N. tabacum W38 plants were grown at 25 °C in long-day conditions (14 h light/10 h dark).

Plasmid construction for transient expression of effectors

genomic DNA was isolated using the Promega Wizard bacterial genomic DNA extraction kit following the manufacturer’s protocol. Allelic variants of RipA1 and RipAA were cloned from South Korean phylotype I isolates as well as GMI1000 and BS048 using Golden Gate [43]. Alleles were PCR-amplified using Takara DNA polymerase (Thermo Fisher Scientific) and cloned into pICH40121 as ~1 kb modules with SmaI/T4 ligase (New England Biolabs). The modules in pICH40121 were assembled into binary vector pICH86988 under control of the constitutive 35S cauliflower mosaic virus promoter and in C-terminal fusion with 6× HA or 3× FLAG epitope. The binary constructs were mobilized into strain AGL1 using electroporation for transient expression experiments in Nicotiana spp.

Hypersensitive response and ion leakage assays

AGL1 carrying effector constructs were grown on LB medium supplied with appropriate antibiotics for 24 h. Cells were harvested and resuspended to final OD600 of either 0.4 (all constructs other than RipH2) or 0.8 (RipH2) in infiltration medium (10 mM MgCl2, 10 mM MES pH 5.6). A blunt syringe was used to infiltrate the abaxial surface of 5-week-old N. benthamiana or N. tabacum leaves. The symptoms of hypersensitive response were observed and photographed at 4 days post-infiltration (p.i.) for N. benthamiana or 6–7 days p.i. for N. tabacum.

Ion leakage assays were performed after infiltration of the same AGL1 strains into N. benthamiana or N. tabacum leaves. Plants were infiltrated with bacterial inoculum (OD600 0.4) and two leaf discs with a diameter of 9 mm each were harvested daily from 0 to 6 days p.i. The leaf discs were placed in 2 ml of distilled water (MilliQ system, Millipore, Bedford, MA, USA) with shaking at 150 r.p.m. for 2 h at room temperature. Ion leakage was measured using a conductivity meter (LAQUAtwin EC-33, Horiba) by loading 60 µl from each sample well. Three or four technical replicates were used for each treatment. Statistical differences between treatments were assessed using Student’s two-tailed t-test with a P value cutoff of 0.05.

Results and discussion

Diversity and distribution of phylotype I isolates across the Republic of Korea

Assemblies of 30 isolates of collected from pepper and tomato in the Republic of Korea between 1999 and 2008 were generated from paired-end Illumina sequencing. Hybrid assemblies were generated with additional Nanopore sequence data for five isolates (Table 1). SPAdes Illumina assemblies were on average 5.7 Mb long, with scaffold counts ranging from 200 (Pe_26) to 306 (To_42). Unicycler hybrid assemblies for five isolates (Pe_1, Pe_3, Pe_27, Pe_39 and Pe_57) were on average 5.85 Mb in length across two to nine scaffolds. Pe_39 and Pe_57 hybrid assemblies are considered complete as they are composed of two to four circularized replicons.

Maximum-likelihood trees were built using non-recombinant core genome alignments generated by readmapping to each Pe_57 replicon as reference. Recombination events identified by ClonalFrameML were removed to obtain non-recombinant trees of South Korean phylotype I strains (Fig. 1). Phylogenies produced using either the chromosome or megaplasmid as reference share nearly identical topologies, with four major clades of isolates and a single representative of a divergent lineage sampled in the eastern mainland (Pe_57). The chromosomal tree is slightly more resolved and has shorter branches leading to clades 1 and 2 than the megaplasmid tree. The two larger clades include strains isolated from both hosts: clade 2 includes isolates from central and northern regions, while clade 4 includes isolates spanning the entire country, including Jeju Island. A sister group to clade 2 only includes pepper isolates (Pe_26, Pe_28, Pe_39 and Pe_1), while a smaller sister group to clade 4 includes a pair of tomato isolates (To_53, To_22). Other than these small clades, there is no strong phylogenetic signal of host of isolation: tomato and pepper isolates are distributed across the phylogeny. Although all but one tomato isolate was sampled from a single northern province and pepper isolates were sampled from seven other provinces, the tomato isolates show no evidence of phylogeographical clustering. Almost every province sampled has isolates from at least two clades, indicating that diverse phylotype I isolates are circulating across many agricultural regions in the Republic of Korea.

Fig. 1.

Diversity and distribution of South Korean phylotype I isolates. (a) Non-recombinant core chromosomal tree and (b) non-recombinant core megaplasmid tree of phylotype I strains isolated from tomato (To) and pepper (Pe) hosts across the Republic of Korea (c). Core alignments were generated by readmapping to separate Pe_57 reference replicons. Phylogenetic trees were built with RAxML v7.2.8 and visualized using FigTree v1.4.3, only displaying branches with a bootstrap support score above 70. Node labels indicate bootstrap support.

Recombination of genes linked with virulence and competition

The cumulative impact of recombination on South Korean phylotype I isolates is enormous: 64.2 % (2 026 835) of the gap-free non-mobile 3 155 880 bp chromosome and 89.0 % (1 511 546) of the gap-free non-mobile 1 698 146 bp megaplasmid alignments are predicted to be affected by recombination in one or multiple strains at some stage in their evolutionary history (Table 2). Recombination occurs at only a slightly reduced frequency relative to mutation on both replicons, but is far more likely to introduce substitutions, with ratios of the effect of recombination/mutation at 4.3 and 5.4 for the chromosome and megaplasmid, respectively. Fewer non-recombinant variant sites were retained for the megaplasmid (614) than the chromosomal (4260) core alignment.

Table 2.

Impact of recombination on South Korean phylotype I isolates

	Alignment length* (bp)	Variable sites	Recombinant sites	R/theta†	1/delta	delta (bp)‡	nu§	r/m||
Chromosome	1 129 045	4260	2 026 835 (64.2 %)	0.805	0.0013	758	0.007	4.30
Megaplasmid	186 600	614	1 511 546 (89.0 %)	0.863	0.0010	1025	0.006	5.40

*Gap-free nonrecombinant alignment length.

†Rate of recombination/mutation.

‡Average length of recombinant regions.

§Average relative divergence of recombinant sequence.

||Relative effect of recombination to mutation.

The location of recombinant events within South Korean phylotype I genomes was assessed to determine whether genes affected by recombination are of relevance to pathogenesis in . Extensive recombination events are present at ancestral nodes in the tree, while the range of predicted recent strain-specific recombination events varies between single recombination events for Pe_15 and Pe_28 to 27 recombination events across both replicons for To_42 and To_63. Since Pe_57 is a unique representative of a separate lineage in the phylogeny, many apparent strain-specific recombination events (118 and 77 on the chromosome and megaplasmid, respectively) are likely to have occurred earlier in the history of this lineage. There are distinct hotspots of recombination on both replicons (Fig. 2, Table S2). Sixty-two chromosomal and 46 megaplasmid genes are affected by two or more recent recombination events, some of which encompass multiple coding sequences (Table S2).

Fig. 2.

Recombination in South Korean phylotype I isolates. Recombination events identified by ClonalFrameML in chromosomal (a) and megaplasmid (b) sequences of South Korean strains. Ancestral recombination events are shown in grey, and recent/strain-specific recombination events are shown in red. Sum of strain-specific events shown at the top of each panel. Genes with two or more strain-specific events have gene IDs listed (hypothetical proteins not listed).

Many genes subject to recent recombination events are linked with virulence and microbial competition. Strikingly, six of these recombination hotspots overlap with six out of seven clusters of haemagglutinin and haemolysin (shlB) genes in Pe_57 (Fig. 2). Some of these clusters appear to be contact-dependent inhibition (CDI) loci, whose products mediate self/non-self-recognition and direct interference competition [44-46]. Contact-dependent inhibition is the product of two-partner (type V) secretion systems that assemble outer-membrane proteins (CdiB) involved in the export and presentation of large, toxic CdiA effectors on the cell surface [47]. CdiA effector protein binding to outer-membrane receptors is followed by the transfer of the C-terminal toxin domain (CdiA-CT) into the susceptible bacterial cell cytoplasm, resulting in cell death. A cognate immunity protein (CdiI) confers resistance to toxin production by neighbouring sister cells by binding and neutralizing specific CdiA-CT toxins; CdiI does not confer immunity to heterologous CdiA-CT toxins. CdiA-CT and CdiI sequences are highly variable both within and between species. Two clusters in Pe_57 (shlB_1 and shlB_4) have cdiA-cdiI-cdiO-cdiB(shlB) arrangements typical of CDI loci in , as well as small ORFs and IS elements between cdiI and cdiB. Comparison of the shlB_4 clusters in isolates for which nanopore assemblies are available (Pe_57, Pe_39, Pe_1, Pe_27 and Pe_3) reveals that Pe_3 diverges significantly from other strains in the C-terminal domain of CdiA and cognate CdiI protein (Fig. S1). Recombination at CDI loci can be a mechanism promoting the maintenance of immunity between otherwise divergent isolates (e.g. Pe_57 vs Pe_39, Pe_1 and Pe_27) or a means of promoting antagonism via rapid diversification in toxin and cognate immunity proteins (e.g. Pe_27 vs Pe_3). The presence of multiple CDI loci in a single strain suggests that both processes may operate simultaneously.

Recombination targets additional genes involved in microbial interactions and virulence. Up to 11 strain-specific recombination events were identified in genes annotated as d-alanine-polyphophoribitol ligase subunit 1 (dltA_4) and tyrocidine synthase tycC_2 genes by Prokka. The closest blastn hit for this region in the phylotype I isolate GMI1000 is to a non-ribosomal peptide synthase locus called ralsomycin (rmyA and rmyB) or ralstonin, which is linked to invasion of fungal hyphae and the formation of chlamydospores [48, 49]. Five chromosomal (ripA1, ripB, ripG6, ripG7, ripTAL) and three megaplasmid (ripH2, ripD, ripAW) effectors show evidence of two or more strain-specific recombination events among the South Korean phylotype I isolates. Flagellar genes (fliI and fliJ), a type IV prepilin-like protein (plpA, annotated by Prokka as pilE1), a type III secretion system translocation protein (hrcU) and component of the export pathway for enterobactin (entS), a high-affinity siderophore, are also affected by multiple strain-specific recombination events in South Korean phylotype I isolates. The detection of recombination within gyrase B and endoglucanase is particularly notable given that these are frequently used to infer relationships between strains.

Core gene tree of the species complex

In order to place the new South Korean strains and T3E dynamics in context of the structure of species complex, a maximum-likelihood phylogeny was inferred from a concatenated alignment of core genes identified in an OrthoFinder pangenome analysis (Fig. 3). The core genome of the 30 South Korean phylotype I genomes and 68 NCBI genomes is represented by 1680 single-copy orthogroups; 17.9 % of the species complex pangenome. The species complex has a large flexible genome, with 39.7 % (3720/9366) of orthogroups present in 15 to 97 genomes and 40.1 % of orthogroups present in the cloud genome (1 to 15 genomes).

Fig. 3.

Phylogeny and type III secreted effector repertoires of the species complex. (a) Maximum-likelihood tree inferred using a non-recombinant core gene alignment of single-copy orthologues. Species and phylotype designation shown above and below branches, respectively. Nodes with bootstrap support values under 70 collapsed. Primary (1°) and secondary (2°) cluster assignment by hierBAPS shown next to strain labels. (b) Type III effector repertoires. Top row displays the normalized Robinson–Fould distance between the effector family gene tree and the core gene tree shown on the left. Lower values (e.g. blue) indicate greater topological similarity between the effector and core gene trees.

The four major phylotypes (and three separate species) are displayed in the core gene tree. The first level of population clustering is consistent with phylotype designations I, IIA, III and IV, while phylotype IIB falls into two separate clusters. Clonal expansion of multiple groups is evident: a phylotype IIB subclade (3.3) typically associated with S. tuberosum infections has been sampled from three continents (South America, Africa and Eurasia). The two larger clades of South Korean phylotype I isolates correspond to subclades 6.2 and 6.4, primarily sampled from solanaceous hosts in PR China, Taiwan and the Republic of Korea (Fig. 3).

The divergent pepper isolate Pe_57 falls into a cluster of East Asian isolates from different hosts: eggplant (EP1, OE1-1); tomato (PSS1308, PSS216, Rs-To2); tobacco (CQPS-1); and squash (RSCM).

Distribution of type III secreted effectors across the RSSC

The distribution of type III secreted effectors across the RSSC was determined to gain an understanding of the extent of variation within Korean isolates and place these in the context of the variation in the RSSC as a whole (Fig. 3b). The average T3E complement varies between phylotypes: phylotype I has on average 71 effector genes per strain, while phylotype II has an average of 61 effector genes per strain. Phylotypes III and IV are each represented by a pair of strains: these have an average of 58 and 64 T3Es per strain, respectively. Many T3E families are distributed broadly across the entire species complex: 64 out of 118 (54 %) of T3Es are present in at least half of the isolates included here. RipAB is present in all strains except YC40-M (Kaempferia sp., PR China) and ripU is potentially pseudogenized in CFBP3059 (Solanum melongena, Burkina Faso), but is otherwise ubiquitous. There are some highly amplified T3E families present across much of the species complex, such as ripA and ripB. Diversification of these families may be linked to host recognition; the broad distribution suggests that they encode critical virulence or host subversion functions.

There are both clade- and phylotype-specific patterns of T3E presence and absence. A set of 14 T3Es (ripA1, ripAG, ripAK, ripBA, ripBE, ripBJ, ripBK, ripBO, ripBP, ripG1, ripS6, ripS8, ripT and Hyp6) are present in 1 or more phylotype I strains, but are absent from phylotype II. This is likewise the case for ripAH, although it is present in two phylotype II strains. Three T3Es (ripJ, ripS5 and ripTAL) are mostly truncated or pseudogenized in phylotype II, while full-length versions are observed in phylotype I. Conversely, seven T3Es (ripBC, ripBG, ripBH, ripBI, ripK, ripS7 and ripV2) and many hypothetical T3Es (Hyp1, Hyp3, Hyp4, Hyp8, Hyp9, Hyp11 and Hyp12) are present in phylotype II but absent in phylotype I.

There are some apparent patterns of recurrent T3E loss across phylotypes I and II. In phylotype II these losses are neither clade- nor host-specific, affecting IPO1609 (S. tuberosum, Netherlands), Po82 (S. tuberosum, Mexico), RS2 (S. tuberosum, Bolivia), UW551 (Geranium sp., Kenya), UW25 (S. lycopersicum, USA), POPS2 (S. lycopersicum, PR China) and RS489 (S. lycopersicum, Brazil). The apparent losses in phylotype I affect tomato strains PSS190 from Taiwan, UTT-25 from India, KACC10709, and To_22 and To_53 from the Republic of Korea. In addition, the Cars-Mep strain isolated from Elettaria cardamom from India shows significant effector losses. These patterns should be interpreted with caution, as some apparent losses (or indeed truncations and pseudogenizations) may be due to variable sequencing and assembly quality.

Despite their ubiquity across the species complex, T3Es are nevertheless subject to recombination. Individual T3E gene trees were compared with the core gene tree to determine whether they exhibit evidence of horizontal transfer. The normalized Robinson–Fould (nRF) distance implemented in ETE3 was used as a measure of topological similarity between individual T3E trees and the core gene tree. Some T3Es, including ripAK, ripBD, ripAH, ripP1, ripC2, ripAG, ripG1 and ripT, have divergent tree topologies (nRF >0.8), indicating that they have been subject to dynamic exchange throughout their evolutionary histories. T3Es with topologies consistent with the core gene tree are typically limited to specific clades (e.g. ripK, ripBO, ripF2 and ripS7, with nRF <0.2). There are, however, three exceptions: ripG5, ripH2 and ripS5 are broadly distributed across the RSSC and are exchanged less frequently, although ripS5 is pseudogenized in many phylotype II strains. Allelic variation is likely to determine fine-scale differences in host range among T3Es with broad distributions across the species complex.

Allelic diversity in T3Es affects host recognition

Allelic diversification of T3Es is likely to play an important role in host adaptation. The link between defence induction and allelic diversity among South Korean isolates was examined for six T3Es: ripA1, ripAA, ripAY, ripE1, ripA5 and ripH2. These families were chosen because they were reported to trigger programmed cell death when overexpressed in Nicotiana spp. [20–22, 24]. Alleles were classified into groups on the basis of shared amino acid sequence identity; patterns of diversity broadly reflected the phylogenetic groupings in the non-recombinant core genome tree for South Korean phylotype I strains (Fig. 4a, b).

Fig. 4.

Effector alleles of RipA1 and RipAA display variable recognition responses in N. benthamiana. Effector alleles were classified into groups (a) based on amino acid sequence alignments. Representatives of each group in RipA1 [(b), top panel] and RipAA [(c), top panel] were tested for their ability to trigger programmed cell death [(b, c), middle panels] and increased ion leakage [(b, c), bottom panels], markers of host recognition. strains expressing allelic variants were infiltrated into N. benthamiana leaves at OD600 0.4 and infiltration sites were photographed 3 days post-infection (p.i.). Asterisks in [(b, c), middle panels) indicate absence of programmed cell death (PCD) for RipA1 allele groups 2 and 4, and RipAA group 2. GFP is a negative control and Pto+AvrPto, GMI1000 RipA1 and BS4048 RipAA refer to positive controls. Ion leakage in [(b, c), bottom panels] was measured at 0, 1, 2 and 3 days p.i. and plotted as the average conductivity of four technical replicates. Error bars represent the standard error of mean (sem).

Representative alleles from each group (Table S3) were transiently expressed in N. benthamiana and N. tabacum plants by agroinfiltration in order to score their ability to elicit programmed cell death and ion leakage, markers of host recognition. Overexpression of RipA1 group 2 and group 4 alleles failed to elicit cell death, while group 1 and group 3 alleles triggered cell death and ion leakage (Fig. 4c, d). One of two substitutions present solely in groups 2 and 4 is likely responsible for the loss of recognition (V1002M); an additional G710D substitution is present in group 4. Host recognition of RipAA alleles varied as well: a single allele (group 2) failed to elicit cell death. A C-terminal truncation preceded by multiple substitutions is likely responsible for the loss of recognition of group 2 RipAA in N. benthamiana (Fig. 4c, d). We found that RipA1 and RipAA protein variants accumulated to significant levels in planta, suggesting that the phenotypes observed are not due to variation in protein stability but rather the loss of host recognition (Fig. S2). Both RipAY alleles present among South Korean phylotype I strains trigger weak cell death and intermediate ion leakage, while all RipE1 and RipA5 alleles caused strong cell death and ion leakage (Fig. S3). In contrast, no RipH2 alleles triggered cell death and ion leakage when overexpressed in N. tabacum (Fig. S3).

Conclusions

We found that South Korean isolates of from tomato and pepper encompass much of the diversity of phylotype I. The location or host from which isolates are obtained have little impact on the population structure of the pathogen in this region. Recent clonal expansions and dissemination between agricultural regions are evident among the South Korean isolates sequenced here. High rates of recombination among the species complex enhance the potential for rapid adaptation and the emergence of novel pathogen variants. Multiple virulence-related genes were affected by recombination events among the South Korean strains. Striking hotspots of recombination overlap with contact-dependent inhibition loci, indicating that microbial competition plays a critical role in evolution, perhaps during stages in its lifecycle when it must compete for nutrients in soil and freshwater environments.

The species complex has an expanded T3E repertoire relative to other bacterial plant pathogens. Many T3Es have evolutionary histories that are inconsistent with the core gene topology, revealing how recombination can reshuffle virulence factor repertoires. We also found evidence of allele-specific recognition of T3Es. Although some effector families may be transferred horizontally less frequently than others, host adaptation and selection imposed by host recognition can result in allelic variation linked to altered host specificity. The integration of pathogen population genomics with molecular plant pathology allows for the identification of conserved pathogen virulence proteins with consistent recognition and elicitation activity, providing suitable targets for breeding crops with durable resistance to phylotype I in the Republic of Korea.

Click here for additional data file.