Introduction

Fusarium oxysporum is a soil-borne fungal pathogen responsible for destructive vascular wilt diseases in a wide range of plants. It ranks within the top ten important fungal pathogens in terms of scientific and economic importance [1]. The best-characterised F. oxysporum pathosystem involves F. oxysporum f. sp. lycopersici (Fol) and tomato. Previous studies of Fol-infected tomato identified a number of fungal proteins within the xylem sap [2]. These secreted in xylem (SIX) effector proteins represent major pathogenicity determinants across different formae speciales of F. oxysporum. Currently, 14 SIX effectors have been identified in Fol consisting of small (less than 300 amino acids in length), secreted, cysteine-rich proteins [3-6]. Most SIX effectors are encoded on the conditionally-dispensable chromosome 14 required for Fol pathogenicity [7]. This dispensable chromosome can be horizontally transferred from Fol to a non-pathogenic strain of F. oxysporum, resulting in a transfer of pathogenicity [4, 7]. To date, all 14 SIX effectors lack sequence identity with proteins of known function, preventing prediction of function based on their amino acid sequence. Several SIX effectors have been shown to be essential for full virulence including SIX1, SIX2, SIX3, SIX5 and SIX6 from Fol [5, 8-11], SIX1 from F. oxysporum f. sp. conglutinans (Focn), which infects cabbage [12], SIX4 from F. oxysporum isolate Fo5176, which infects Arabidopsis [13], and SIX1 and SIX8 from F. oxysporum f. sp. cubense, which infects banana [14, 15]. Fol SIX3 (Avr2) and SIX5 are adjacent, divergently-transcribed genes with a common promoter, and SIX5 has been shown to interact with SIX3 to promote virulence by enabling symplastic movement of SIX3 via plasmodesmata [16]. Focn SIX8 and PSE1 (pair with SIX8 1) are also a divergently-transcribed effector gene pair that function together to suppress phytoalexin production and plant immunity in Arabidopsis [17]. In Fol, SIX8 forms a similar gene pair with PSL1 (PSE1-like 1) [17]. Despite their roles in fungal pathogenicity, the virulence functions of most SIX effectors remain unknown.

To combat pathogen attack, plants possess resistance genes that encode immunity receptors capable of recognising specific effectors leading to disease resistance. Four resistance genes, introgressed into tomato from related wild species, have been cloned. I and I-7 encode transmembrane receptor proteins containing extracellular leucine-rich repeat (LRR) domains and short cytoplasmic domains (LRR-RPs) [18, 19]. I-2 encodes a cytoplasmic receptor containing nucleotide binding (NB) and C-terminal LRR domains [20], while I-3 encodes a transmembrane protein with an extracellular S-receptor-like domain and cytoplasmic serine/threonine kinase domain (SRLK) [21]. Fol Avr1 (SIX4), Avr2 (SIX3) and Avr3 (SIX1) are recognised by tomato immunity receptors, I, I-2 and I-3, respectively, leading to effector-triggered immunity (ETI) and disease resistance [5, 22, 23].

By understanding the function of F. oxysporum effector proteins, and how specific effectors are detected by immunity receptors, we (and others) hope to develop novel disease management strategies targeting vascular wilt diseases. Protein structure studies of effectors provide one avenue to assist this pursuit. Currently, Avr2 represents the only SIX effector whose protein structure has been determined [24]. Interestingly, the β-sandwich fold of Avr2 revealed that this effector shares structural homology to ToxA from Pyrenophora tritici-repentis and AvrL567 from Melampsora lini [25, 26], despite a lack of sequence identity. The observation of structural classes for effectors without identifiable domains or homologies to proteins of known function has been demonstrated experimentally for four effector structural families, including the so-called MAX (Magnaporthe oryzae Avr effectors and ToxB from P. tritici-repentis) [27], RALPH (RNAse-Like Proteins associated with Haustoria) [28], LARS (Leptosphaeria Avirulence and Supressing) [29] and ToxA-like families [24-26].

Combining experimental and computational approaches, we present the structural repertoire of sequence unrelated effectors utilised by Fol during infection of tomato, including the classification of a new effector family, the FOLD (Fol dual-domain) effectors. We show using structural comparisons that FOLD effectors are widely distributed in phytopathogenic fungi as well as symbionts. Further, we define the domains and residue that mediate the recognition of the FOLD effector, Avr1, by its corresponding immunity receptor.

Results

The structures of Avr1 and Avr3 adopt a similar dual-domain fold

Avr1 and Avr3 are cysteine-rich effectors that belong to the K2PP (Kex2-processed pro-domain) effector class [30, 31]. To help understand their function, and recognition by I and I-3, we sought to solve their structures using x-ray crystallography. Using our optimised protein production strategy [32], we produced Avr1 (Avr1^18-242) and Avr3 (Avr3^22-284) in Escherichia coli for crystallisation studies (S1A and S1B Fig). Crystals were obtained for Avr3^22-284 (after referred to as Avr3) (S1B Fig), however, Avr1^18-242 failed to crystallise. Previously, we demonstrated that pro-domain removal from the K2PP effector SnTox3 was required to obtain protein crystals [30] and predicted this may also be important for Avr1. Treatment of Avr1^18-242 with Kex2 in vitro resulted in a predominant Avr1 band of ∼20 kDa consistent with a mature Avr1^59-242 protein, however, lower molecular weight bands were also observed suggesting in vitro Kex2 cleavage at additional sites [30]. To address this, Avr1 was engineered with an internal thrombin cleavage site (replacing the Kex2 site) to produce a single Avr1^59-242 product after thrombin cleavage (after referred to as Avr1). This protein was subsequently used for crystallisation studies resulting in rectangular plate-like crystals (S1A Fig).

The crystal structures of Avr1 and Avr3 were solved using a bromide-ion-based single-wavelength anomalous diffraction (SAD) approach (S1 Table), and subsequently were refined using a native dataset to a resolution of 1.65 Å and 1.68 Å, respectively (Fig 1A and 1B). Despite sharing low amino-acid sequence identity (19.5%), Avr1 and Avr3 adopt a structurally similar dual-domain protein fold. Interpretable, continuous electron density was observed from residue 96 in Avr3 and some regions of the intact pro-domain could be interpreted in the electron density (residues 26-49) (S2A Fig). We also identified regions of the pro-domain (residues 23-45) of Avr1 in the electron density, despite thrombin cleavage of the pro-domain prior to crystallisation (S1A Fig). This indicates that an association between respective Avr and pro-domain was maintained post cleavage in vitro (S2B Fig). The importance of this association, if any, remains unclear, but for simplicity, the pro-domains were excluded from further analysis.

Fig 1.

The Avr1 and Avr3 N-terminal domain (N-domain), consisting of an N-terminal α-helix followed by five β-strands, and C-terminal domain (C-domain), consisting of a β-sandwich architecture, involving seven or eight β-strands are very similar with a root-mean-square deviation (RMSD) of 2.1 Å and 2.8 Å, respectively (superposition performed using DALI server [33]) (Fig 1). While the individual domains are very similar, superposition of the dual-domain structures returns an RMSD of ∼3.4 Å. The larger difference is due to a rotation between the N-and C-domains (Fig 1E). The structures of Avr1 and Avr3, when compared with the solved structures of other fungal effectors, demonstrate that they adopt a unique two-domain fold and represent the founding members of a new structural class of fungal effectors we have designated the FOLD (Fol dual-domain) effectors.

SIX6 and SIX13 belong to the FOLD effector family

We were interested to determine if other SIX effectors belonged to the FOLD effector family. One conserved sequence feature observed in Avr1 and Avr3 was the spacing of the six cysteines within the N-domain. We analysed the cysteine spacing of the other SIX effectors and found that SIX6 and SIX13 contained a cysteine profile like Avr1 and Avr3 (Fig 2A), suggesting they may be FOLD effectors. With the recent advances in ab initio structural prediction by Google DeepMind’s AlphaFold2 [34] we predicted the structures of the SIX effectors to determine if, as suggested by our sequence analysis, other SIX effectors are FOLD effector family members.

Fig 2.

As an initial step we benchmarked the AlphaFold2 predicted models of Avr1 and Avr3 (downstream of the Kex2 cleavage site (Avr1^59-242 and Avr3^96-284) against our experimentally determined structures (S3 Fig). The AlphaFold2 model of Avr1 returned a low average per-residue confidence score (pLDDT = 55%) and the RMSD was 6.9 Å when model and structure were compared, however, the dual domain architecture was correctly predicted with a Z-score of 11.3 identified using a DALI pair-wise structural comparison (S3A Fig and S3E). The AlphaFold2 model of Avr3 returned a high average pLDDT score (92%) and superimposed well to the solved structure (S3B Fig), despite a slight skew between the orientation of the individual domains (RMSD = 3.6 Å overall; 1.1 Å for the N-domain; 0.8 Å for the C-domain). This demonstrated that accurate FOLD effector prediction was possible using AlphaFold2.

We subsequently generated SIX6 and SIX13 models, downstream of the predicted Kex2 cleavage site (SIX6^58-225, SIX13^78-293), using AlphaFold2 and obtained high average confidence scored models supporting their inclusion in the FOLD family (S4 Fig). To validate this experimentally, we produced SIX6^58-225 and SIX13^22-293 (after referred to as SIX6 and SIX13) as described for Avr1/Avr3 and obtained crystals for both proteins (S1 Fig). While the SIX13 crystals diffracted poorly, the SIX6 crystals diffracted x-rays to ∼1.9 Å and we solved the structure of SIX6 using the AlphaFold2 generated model as a template for molecular replacement (Fig 2B, S1 Table), confirming its inclusion as a member of the FOLD family. Despite lacking an N-terminal helix, the N-domain contains five β-strands held together by three disulfide bonds with an arrangement, identical to Avr1 and Avr3. The C-domain is an eight stranded β-sandwich that is stabilised by a single disulfide bond (unique to SIX6 compared to Avr1 and Avr3) connecting the β7 and β12 strands. Like Avr1, we identified regions of the pro-domain within the SIX6 structure (residues 29-46), despite cleavage of the pro-domain prior to crystallisation (S2C Fig), but only within one molecule in the asymmetric unit (S2D Fig). For structural analysis, we used the structured region of Chain A of SIX6 (Fig 2B).

FOLD effectors are distributed across multiple fungal genera

Despite structural similarities, the FOLD effectors are divergent in their amino acid sequences, sharing 15.5 – 22.5% sequence identities between all members (Fig 2A). Homologues of FOLD effectors are dispersed across multiple formae speciales of F. oxysporum (S5 Fig) [6, 8, 35-38]. We were interested to understand the distribution of FOLD effectors in fungi. Previous structural-based searches performed on effector candidates from Venturia inaequalis using Avr1 and Avr3 as templates (which we provided to the authors) found three candidates predicted to be FOLD effectors [39]. Here, we utilised our experimentally determined structures (Avr1, Avr3 and SIX6) to search for other fungal FOLD effectors within the AlphaFold2 protein structure database [40] (https://alphafold.ebi.ac.uk/) using the Foldseek webserver [41]. This analysis identified 124 putative FOLD protein family members across three Divisions of Fungi (Ascomycota, Basidiomycota, and Glomeromycota) (Fig. 2C). Over half of these were found in Ascomycota fungi (73), with expanded families in species of Fusarium and Colletotrichum (Fig 2C, S2 Table). Expanded families of FOLD proteins were also observed in the division Glomeromycota that form arbuscular mycorrhiza in plant roots, while two putative FOLD effectors were also predicted in the ectomycorrhizal fungus Piloderma olivaceum (Division Basidiomycota), which forms mutualistic associations with conifer and hardwood species [41]. Structural superposition of members from the three Divisions confirms the structural similarities between the N and C domains and highlights that the major differences identified are the orientation of the domains relative to each other (Fig. 2D), consistent with our experimental data for Avr1, Avr3 and SIX6.

Distinct structural families exist among the other SIX effectors

With the successful utilisation of AlphaFold2 as a model for molecular replacement (SIX6 structure), and structural similarity searches for FOLD effectors, we decided to perform structural comparisons with the remaining SIX effectors. AlphaFold2 modelling of the effectors was conducted on sequences with the signal peptide and putative pro-domain (if present) removed (S6 Fig). The models and experimentally determined SIX effector structures (Avr1, Avr2, Avr3 and SIX6) were compared using the DALI server [33] and a Z-score with a cutoff of >2 was used to indicate structure similarity.

The observed structural similarity between the FOLD effectors was high, with Z-scores above 8 for all comparisons (Fig 3A). Avr2, a member of the ToxA-like effector family, exhibited structural similarity with the SIX7 and SIX8 models (Z-scores > 5) (Fig 3A). Analysis of the models and topology show that SIX7 and SIX8 both consist of a β-sandwich fold, strongly indicating their inclusion of within the ToxA-like structural family (Fig 3C, S7 Fig).

Fig 3.

Beyond these described structural families, the Z-scores indicated that two additional, but not yet characterised, structural families exist within the SIX effectors. Here, we define these as structural family 3 and 4, consisting of SIX9 and SIX11, and SIX5 and SIX14, respectively (Fig 3D, E). The models of SIX9 and SIX11 both consist of five β-strands and either two or three α-helices (Fig 3D, S8 Fig), despite sharing only 14% sequence identity. To further our understanding of the putative function of this family we did a structural search against the protein databank (PDB) and found that both structures share structural similarity to various RNA binding proteins (Z-scores > 2.5), such as the RNA recognition motif (RRM) fold of the Musashi-1 RNA-binding domain (PDB code: 5X3Z) [42] with a Z-score of 2.6 and 4.5 for SIX9 and SIX11, respectively.

SIX5 and SIX14 also share limited sequence identity (23%) but the structural predictions show a similar secondary-structure topology consisting of two α-helices and four to six β-strands (Fig 3E, S8 Fig). We compared the models of SIX5 and SIX14 against the PDB using DALI and identified structural similarity toward the Ustilago maydis and Zymoseptoria tritici KP6 effector (PDB codes: 4GVB and 6QPK) [43], suggesting SIX5 and SIX14 belong to the KP6-like structural family (S7 Fig). Collectively, this analysis demonstrates that 11 of the 14 SIX effectors, group into 4 different structural families.

Structural modelling and comparison of an expanded set of Fol effectors

The SIX effectors are only a subset of effectors utilised by Fol during infection of tomato. Recently, the Fol genome was re-sequenced [44] and reannotated in combination with RNAseq data from Fol- infected tomato plants [45]. A total of 26 genes encoding novel effector candidates were identified that were consistently upregulated during Fol infection [45], which were not previously predicted or predicted incorrectly in the original genome annotation [4]. Of these, 14 genes encoded proteins with no recognised domains or motifs based on their amino acid sequences. We generated structural models using AlphaFold2 of these 14 (S3 Table, S6 Fig) and structurally aligning them using DALI against SIX effector representatives from each family to assess if they fell into any of the established families (Fig 3B). We found the predicted structure of FOXGR_015533 adopts a nine β-stranded sandwich and is likely a member of the ToxA-like class (Fig 3C). PSL1 [17] and FOXGR_015322, here designated PSL2, are sequence related effectors (∼85% sequence identity) and show a conserved structure (Fig 3E). Both have Z-scores of >2 against Family 4 and are likely members of this family.

Based on this analysis we also suggest an additional structural family. FOXG_18699 and FOXGR_015522 are structurally related (Z-score of 2.2) with a sequence identity of ∼29%. While FOXGR_015522 does share some resemblance to Family 4, based on manual alignment (Fig 3F) and domain topology analysis (S8 Fig) these effectors appear to belong to an independent structural family, designated Family 5. Collectively, these data demonstrate that Fol utilises multiple structurally related, sequence diverse, effectors during infection of tomato.

Interaction between effector pairs from two structural families

In Fol, Avr2 and SIX5, and SIX8 and PSL1 form a similar head-to-head relationship in the genome with shared promoters and are divergently-transcribed (Fig 4A) [16, 17]. Previously, studies concerning Avr2 and SIX5 have demonstrated that the proteins function together and interact directly via yeast-two-hybrid analysis [9]. Homologues of SIX8 and PSL1 from Focn (SIX8 and PSE1) are also functionally dependent on each other, however an interaction could not be established in yeast [17]. Here we demonstrate that both protein pairs contain a ToxA-like family member (Avr2, SIX8) and a structural family 4 member (SIX5, PSL1). Considering the predicted structural similarities, we were interested in testing whether Fol SIX8 and PSL1 interact.

Fig 4.

We heterologously produced Fol SIX8^50-141 and PSL1^18-111 (S1E-F Fig) (after referred to as SIX8 and PSL1) and co-incubated the proteins before analysing by size exclusion chromatography (SEC) (Fig 4B). The elution profile of PSL1 alone showed a major peak (∼12.25 mL) at a volume consistent with a dimeric form of the protein, while SIX8 showed a major peak (∼15 mL) consistent with a monomer (Fig 4B). Strikingly, when incubated together the major protein peaks migrate to ∼12.8 mL. SDS-PAGE analysis confirmed that presence of PSL1 and SIX8, indicating that the migration of both proteins on SEC is altered after incubation (Fig 4B). These data are consistent with PSL1 and SIX8 forming a heterodimer.

To understand the structural basis of the interaction, we attempted to solve the structure of the complex, but we were unable to obtain crystals. We subsequently utilised AlphaFold2-Multimer [46] through ColabFold [47], to model the interaction. Manual inspection of the top 5 models (S10A Fig, top model shown Fig. 4C) demonstrated that the thiol side chain of a free cysteine in PSL1 (Cys 37) and SIX8 (Cys 58) co-localised in the dimer interface, suggesting that an inter-disulfide bond may mediate the interaction. To test this, we performed intact mass spectrometry of SIX8 and PSL1 (alone and post incubation) under non-reduced and reducing conditions. The mass observed from the incubated SIX8 and PSL1 non-reduced sample contained a predominant species consistent with the combined molecular weight of SIX8 and PSL1 (20777 Da) (Fig 4D, S9G-H Fig). SIX8 and PSL1 failed to form a heterodimer with an unrelated protein containing a free cysteine, suggesting specificity in the interaction (S9I-L Fig). Collectively, these data demonstrated that the SIX8-PSL1 heterodimer is mediated via a disulfide bond.

To confirm the involvement of the predicted residues involved, interaction with cysteine mutants of PSL1 and SIX8 (PSL1_C37S^18-111 and SIX8_C58S^50-141, after referred to as PSL1_C37S and SIX8_C37S) were analysed (Fig 4E). When PSL1_C37S was incubated with SIX8_C37S or SIX8 alone, the heterodimer was not resolved via SEC (Fig 4D, S10B Fig). This was further confirmed using mass spectrometry (Fig 4C). We crystallised and solved the structure of SIX8_C58S at 1.28 Å (S1E Fig and S10C Fig) which confirms its inclusion within the ToxA-like structural family (S10D Fig).

The molecular basis of Avr1 recognition by the I receptor

The structural identification of the FOLD effector family provides an opportunity to understand their recognition by cognate immunity receptors. Here, we focussed on Avr1 (SIX4), which is recognised by the I immunity receptor leading to ETI and disease resistance [19]. Previous studies have shown co-expression of the I gene from the M82 tomato cultivar (I^M82) with Avr1 in Nicotiana benthamiana leads to a cell death response, a proxy for ETI [19]. Conversely, co-expression with the allelic variant (i^Moneymaker) from the susceptible cultivar Moneymaker does not lead to cell death as the receptor cannot recognise Avr1 [19] (Fig 5B). Here we sought to further define the recognition between Avr1 and I utilising the N. benthamiana system.

Fig 5.

To facilitate this, we identified homologues of Avr1 that possess natural residue variation. FonSIX4, a homologue of Avr1 from the watermelon pathogen, F. oxysporum f. sp. niveum (Fon) shares 79% identity with Avr1 (Fig 5A). Using the N. benthamiana assay we show FonSIX4 is recognised by I receptors from both cultivars (I^M82 and i^Moneymaker) (Fig 5B). FonSIX4 and Avr1 differ by 34 residues distributed across both N-and C-domains of the protein (Fig 5A). To narrow down the regions involved in recognition we designed chimeric variants by swapping the N-and C-domains (Avr1^NFonSIX4^C and FonSIX4^NAvr1^C) (Fig 5A, C). When these were co-expressed with i^Moneymaker the cell death response, quantified using ion leakage assays (Fig 5D-E) and visual inspection (S11A Fig), suggest the C-domain of FonSIX4 is recognised by i^Moneymaker. We separated Avr1 and FonSIX4 proteins into their N-or C-domains and co-expressed with I^M82 or i^Moneymaker. Quantification using ion leakage assays demonstrate that the C-domains of Avr1 and FonSIX4 cause cell death when expressed with I^M82 and I^{M82/Moneymaker}, respectively. These data confirm the C-domain is sufficient for I receptor recognition (Fig 5D-E, S12 Fig, see S11 Fig for N. benthamiana leaf infiltration and protein accumulation data).

To understand how Avr1 can escape i^Moneymaker recognition, we focussed on surface exposed variant residues (underlined) mapping to four regions within the C-domain (Fig 5A and 5C). Four reciprocal swap mutants between Avr1 and FonSIX4 (Avr1^ADVKT, Avr1^IDH, Avr1^NGQAR, Avr1^EEEYGIN) were co-expressed with i^Moneymaker to identify the residues required for FonSIX4 recognition. Avr1^EEEYGIN showed consistent ion leakage and cell death similar to FonSIX4 (Fig 5F-G), whereas ion leakage quantification for the other three mutants (Avr1^ADVKT, Avr1^IDH, Avr1^NGQAR) was statistically similar to the non-recognised Avr1 (Fig 5G). The reciprocal mutations in FonSIX4 (FonSIX4^KEVYHID) significantly reduced ion leakage and cell death response when co-expressed with i^Moneymaker compared to FonSIX4 (S11D-E Fig, see S11G Fig for protein accumulation data). Collectively, these data show that the C-domain in Avr1 is recognised by I^M82, and surface exposed residues in the C-domain allow Avr1 to escape recognition by i^Moneymaker.

Discussion

Pathogenic fungi are in a continuous arms race with their plant hosts. To aid virulence, but avoid detection, effectors evolve rapidly causing significant diversity at the amino acid sequence level [48]. An emerging theme in fungal effector biology is the classification of effectors into families based on structural similarity [49]. Here, we demonstrate that despite their sequence diversity, the Fol SIX effectors can be classified into a reduced set of structural families. This observation has implications for functional studies of SIX effectors, and ultimately our understanding of the infection strategies used by F. oxysporum.

Expanding the structural classes in fungal effectors

To date, five fungal effector families have been defined based on experimentally-determined structural homology, including the MAX [27], RALPH [28, 50, 51], ToxA-like [24-26], LARS [29, 52] and FOLD effectors, defined here. Effectors that fall within many of these structural families are shared across distantly related fungal species. The ToxA-like family includes effectors from fungi that group to both divisions of higher-fungi (Basidiomycota and Ascomycota fungi) [24-26]. The MAX effector family were originally defined as AVR effectors from M. oryzae and ToxB from P. tritici-repentis [27] but pattern-based sequence searches suggest they are widely distributed amongst the Dothideomycetes and Sordariomycetes [27, 53]. Similarly, LARS effectors, defined in Leptosphaeria maculans and Fulvia fulva, have structural homologues predicted in at least 13 different fungal species [29]. Based on sequence homologues alone, FOLD effectors are well dispersed in fungi with homologues amongst the Sordariomycetes including many formae speciales of F. oxysporum, Colletotrichum and Ustilaginoidea. Based on structural comparison of the AlphaFold2 structural database we show that is extended to fungi in three Divisions, including plant pathogens and symbionts. This was supported by a recent study modelling the secretomes of arbuscular mycorrhizal fungi which found enlarged and diversified gene families encoded proteins predicted to share the FOLD effector structure [54]. The exclusive presence of FOLD effectors in plant-colonising fungi may suggest they facilitate plant colonisation in pathogenic and symbiotic fungi [54].

Effector structure prediction

Experimentally determining the structures of fungal effectors is not a trivial undertaking. From challenges associated with effector protein production through to hurdles related to structure solution (such as experimental phasing), the research time required to determine an effector structure experimentally ranges from months to many years (sometimes never). Not surprisingly, any reliable structural modelling methods are welcomed by researchers interested in effector biology. To this end, several recent studies have used effector structure prediction to expand our understanding of plant-microbe interactions [55, 56].

Work by Bauer and colleagues, prior to the release of AlphaFold2, used structural modelling to show that numerous recognised Avr effectors from the barley powdery mildew-causing fungal pathogen Blumeria graminis (Bgh) are members of the RALPH effectors class [55]. Seong and Krasileva used similar structural modelling approaches to predict the folds of ∼70% of the Magnaporthe oryzae secretome [56]. In doing so, they suggested an expansion in the number of MAX effectors and identified numerous sequence-unrelated groups of structural homologues (putative structural classes) within M. oryzae. Making use of AlphaFold2, Yan and colleagues show that structurally conserved effectors, including the MAX effector family, from M. oryzae are temporally co-expressed during the infection process [57]. In the largest comparison study to date, Seong and Krasileva carried out a large comparative structural genomics study of fungal effectors utilising AlphaFold2 [58]. Their findings support the hypothesis that the structurally conserved effector families are the result of divergent evolution and support previous finding that the structural landscape of effectors is more limited than what is suggested by sequence-diversification.

Here, we were in a unique position to apply and benchmark AlphaFold2 against experimentally determined structures for Fol effector prediction. We subsequently used AlphaFold2 to demonstrate that, within the repertoire of effectors we tested, up to five sequence-unrelated structural families are secreted during Fol infection. There are numerous caveats in relying solely on AlphaFold2 to generate structural models of effectors. The accuracy of models generated by AlphaFold2 can decline in cases with low numbers of homologues (∼30 sequences in the multiple sequence alignment) [34]. This may help explain the low confidence prediction for SIX4 (Avr1) (S4A Fig), which is only distributed in a few ff. spp. of F. oxysporum. This poses a potential issue for predicting the structures of fungal effectors that lack homologues. In our hands, we have had mixed results when comparing several unpublished effector structures experimentally determined in our lab to AlphaFold2 models. In some instances, the models are wrong, for example AvrSr50 [59], however, in these cases the AlphaFold2 predictions reported low confidence scores, an important criterion for assessment of model reliability. Despite this, AlphaFold2 models were critical in solving the structure of SIX6 and SIX8, as templates for molecular replacement. This negated the need to derivatise our crystals, a process that we had struggled with for SIX6 crystals, significantly reducing the time and research effort to determine the experimental structures.

Structural classes: A starting point for functional characterisation

Given their lack of sequence identity to proteins of known function or conserved motifs, structural determination of effectors is often pursued to provide functional insight and understanding of residues involved in recognition. The existence of structural families of effectors raises the question of whether links can now be made concerning their function based on structural similarities. Unfortunately, the FOLD effectors share little overall structural similarity with known structures in the PDB outside of the similarity with each other. At a domain level, the N-domain of the FOLD effector Avr1 has some structural similarities with cystatin cysteine protease inhibitors (PDB code: 4N6V, PDB code: 5ZC1) [60, 61], and the C-domain with tumour necrosis factors (PDB code: 6X83) [62] and carbohydrate-binding lectins (PDB code: 2WQ4) [63]. Relatively weak hits were observed for Avr3/Six6.

Interestingly, the predicted models for SIX9 and SIX11 within Family 3 have structural homology with RNA-binding proteins (PDB code: 3NS6, PDB code: 5X3Z) [42, 64], unrelated to RALPH effectors. Despite this structural homology, close inspection of these models suggests RNA binding is unlikely, as in both models the putative RNA binding surface is disrupted by a disulfide bond.

The putative family 4 effectors (SIX5, SIX14, PSL1 and PSL2) have structural homology with KP6 effectors and heavy metal associated (HMA) domains. Metal binding within HMA domains is facilitated by conserved cysteine residues [65], however, their absence in the family 4 effectors suggests they are unlikely to have this activity.

The putative family 5 effectors (FOXGR_015522 and FOXG_18699) have structural homology with different proteins within the PDB. FOXGR_015522 is structurally similar to plant defensins (PDB code: 6MRY, PDB code: 7JN6) [66, 67] and K⁺ channel-blocking scorpion toxins (PDB code: 1J5J, PDB code: 2AXK) [68, 69]. FOXG_18699 has structural homology with the C-terminal domain of bacterial arginine repressors (PDB code: 1XXB, PDB code: 3CAG) [70, 71].

A structural explanation for functional effector pairs

One interesting outcome of this study is a link between structural families and co-operative interactions between effectors. The ToxA-like effectors, Avr2 and SIX8 are known to form functional effector pairs with SIX5 and PSE1 (PSL1-homolouge), respectively [9, 17]. According to our modelling work, both SIX5 and PSL1 are members of structural family 4. Avr2 and SIX5 are adjacent divergently-transcribed genes on Fol chromosome 14 and the protein products have been shown to physically interact [9]. Likewise, SIX8 and PSL1 are adjacent divergently-transcribed genes in the Fol genome and we demonstrate here a physical interaction between the proteins. The AlphaFold2-multimer models of the SIX8 and PSL1 heterodimer, drew our attention to the inter-disulfide bond between SIX8 and PSL1 required for the interaction, which we confirmed experimentally. While these residues are conserved in Focn SIX8 and PSE1, the Avr2 structure and SIX5 model lack free cysteine residues, suggesting a different mode of interaction.

Interestingly, two other SIX genes also form a divergently-transcribed gene pair on Fol chromosome 14. SIX7 (ToxA-like family) and SIX12 possess start codons 2,319 base-pairs apart and potentially share a common promoter. While SIX12 did not group with any structural families, the AlphaFold2 model had a very low prediction confidence (35.5%). On closer inspection of the sequence, we observed that the cysteine spacing in SIX12 closely resembles other family 4 members (S13 Fig), which suggests that SIX12 may also be a family 4 member. We therefore speculate that SIX7 and SIX12 may function together, as described for the Avr2/SIX5 and SIX8/PSL1 pairs.

Are experimentally derived effector structures still worth the effort?

The potential of machine-learning structural-prediction programs, such as AlphaFold2, heralds an exciting era, especially for a field that has long suffered from a lack of prediction power based on effector sequences. A question now emerges; when prediction model confidence is high, should we bother solving structures experimentally? The answer to such a question will always depend on what the structure is being used for. Ultimately, structural models, whether experimentally or computationally derived, represent information to base and/or develop a hypothesis to subsequently test. Here we demonstrate the power of structure prediction in combination with experimentation, both for validating models and understanding protein:protein interaction interfaces. One interesting observation we made was that while the AphaFold2-multimer models of the SIX8 and PSL1 heterodimer were sufficient to highlight the cysteine residues required for mediating the interaction, the models and interaction interfaces differed significantly (S10A Fig). When the modelling was repeated with the SIX8^C58S experimentally derived structure included as a template, the interaction models and heterodimer interface were of higher quality and essentially identical (S10E Fig). This observation can be retrospectively reconciled. The region of SIX8 involved in the interaction with PSL1 was modelled incorrectly by AlphaFold2 when compared to the structure (S10D Fig). Collectively, these data highlight that some models are good enough, but others maybe better.

Effector structural classes and understanding receptor recognition

Understanding the structural basis of plant immunity receptor-effector interactions represent a key step towards engineering plant immunity receptors with novel specificities. Recent structures of nucleotide-binding domain leucine-rich repeat (NLR) proteins reveal exquisite details of these direct interactions [72-74]. The FOLD effectors, Avr1 and Avr3, are recognised by different classes of immunity receptors; I, an LRR-RP [19] and I-3, a SRLK [21]. While the mode of recognition has not yet been described for Avr3, we demonstrate here that Avr1 is recognised at the C-domain (Fig 5). This is significant because it demonstrates that different immunity receptor classes can recognise structural homologues. It might also help explain the function of Avr1 during Fol infection. When Houterman and colleagues identified Avr1, they demonstrated that it could suppress plant immunity conferred by the I-2 and I-3 receptors [22]. Considering our structural understanding of these FOLD effectors, it is plausible that Avr1 achieves suppression of I-3-mediated immunity by preventing Avr3 recognition through competitive inhibition. The LARS effectors represent another example of effectors that can activate and suppress resistance-gene-mediated immunity. AvrLm4-7 can prevent recognition of AvrLm3 and AvrLm9 (all LARS structural homologues [29]), by their cognate Rlm receptors [75, 76]. Rlm4, Rlm7 and Rlm9 all encode for wall-associated kinases [77, 78], but the identify of Rlm3 remain unknown. These studies demonstrate that members of at least two different structural effector families can suppress immunity triggered by structurally homologous effectors.

Collectively, the results presented here will aid future studies to understand the molecular basis of F. oxysporum effector function and recognition, and by extension, the design and engineering of immunity receptors with novel recognition specificities to help protect plants against Fusarium wilt disease.

Materials and methods

Vectors and gene constructs

SIX6, Avr1Thrombin, SIX6-TEV, SIX8Thrombin, SIX8_C58SThrombin, PSL1, PSL1_C37S and SIX13 coding sequences (without their signal peptides as determined by SignalP-5.0) were codon optimised for expression in E. coli and synthesised with Golden-Gate compatible overhangs by Integrated DNA Technologies (IDT, Coralville, USA) (S4 Table). The Kex2 cleavage motif of Avr1 and SIX8 were replaced with a thrombin cleavage motif, and TEV protease cleavage motif for SIX6 for pro-domain processing. Avr1 and Avr3 coding sequences were PCR amplified using Fol cDNA as a template with primers containing Golden-Gate compatible overhangs. All genes for E. coli expression were cloned into a modified, Golden-Gate-compatible, pOPIN expression vector [79]. The final expression constructs contained N-terminal 6xHis-GB1-tags followed by 3C protease recognition sites. The Golden-Gate digestion, ligation reactions and PCR were carried out as described by Iverson, Haddock [80]. Avr1 and FonSIX4 mutant sequences without the signal peptide were synthesised with compatible overhangs by IDT (S4 Table) and cloned into the pSL vector containing the Nicotiana tabacum PR1 signal peptide using the In-fusion cloning kit (Takara Bio USA Inc., San Jose, USA) to allow efficient secretion in N. benthamiana via Agrobacterium-mediated expression. For tagged constructs, Avr1 and FonSIX4 mutant sequences and 3xHA tag were amplified with PCR and assembled using In-fusion cloning into the pSL vector. All of the primers were synthesised by IDT (S5 Table). All constructs were verified by sequencing.

Protein expression and purification

Sequence-verified constructs were co-expressed with CyDisCo in SHuffle T7 Express C3029 (New England Biolabs (NEB), Ipswich, USA) and purified as previously described [32]. For Avr3, the buffers used after fusion tag cleavage were altered slightly to increase protein stability and a second IMAC step was excluded after the cleavage of the N-terminal fusion tag. During the cleavage step, the protein was dialysed into a buffer containing 10 mM MES pH 5.5 and 300 mM NaCl. The size-exclusion chromatography (SEC) HiLoad 16/600 Superdex 75 pg column (Cytiva) was equilibrated with a buffer containing 10 mM MES pH 5.5 and 150 mM NaCl.

For biochemical and crystallisation studies, Avr1 and SIX8 with an internal thrombin cleavage site, and SIX6 with an internal TEV protease cleavage site for pro-domain removal were processed with 2 to 4 units of thrombin from bovine plasma (600-2,000 NIH units/mg protein) (Sigma-Aldrich Inc., St. Louis, USA) per mg of protein at 4°C or TEV protease (produced in-house) until fully cleaved. Fully cleaved proteins were purified further by SEC using a HiLoad 16/600 or HiLoad 26/600 Superdex 75 pg column (Cytiva) equilibrated with a buffer containing 10 mM HEPES pH 7.5 or 8.0 and 150 mM NaCl. Proteins were concentrated using a 10 or 3 kDa molecular weight cut-off Amicon centrifugal concentrator (MilliporeSigma, Burlington, USA), snap-frozen in liquid nitrogen and stored at -80°C for future use.

Intact mass spectrometry

Proteins were adjusted to a final concentration of 6 µM in 0.1% (v/v) formic acid (FA) for HPLC-MS analysis for untreated samples. For reduced samples, DTT was added to the protein to a final concentration of 10 mM. Proteins were incubated at 60°C for 30 minutes and adjusted to 6 µM in 0.1% (v/v) FA. Intact mass spectrometry on all proteins was carried out as described previously [32]. Data were analysed using the Free Style v.1.4 (Thermo Fisher Scientific) protein reconstruct tool across a mass range of m/z 500 – 2000 and compared against the theoretical (sequence based) monoisotopic mass.

Circular dichroism (CD) spectroscopy

The CD spectra of purified effectors of interest were recorded on a Chirascan spectrometer (Applied Photophysics Ltd., UK) at 20°C. Samples were diluted to 10 µM in a 20 mM sodium phosphate buffer at pH 8.0. Measurements were taken at 1 nm wavelength increments from 190 nm to 260 nm. A cell with a pathlength of 1 mm, a bandwidth of 0.5 nm and response time of 4 s were used, with 3 accumulations. The data were averaged and corrected for buffer baseline contribution, and visualised using the webserver CAPITO tool with data smoothing [81]. CD analysis was performed on all purified proteins (S14 Fig).

Crystallisation, diffraction data collection and crystal structure determination

Initial crystallisation screening was performed for Avr3^22-284, Avr1^18-242, Avr1^59-242, SIX8^50-141, PSL1^18-111, SIX6^17-225, SIX6^58-225, SIX8_C58S^19-141, SIX8_C58S^50-141, PSL1_C37S^18-111, SIX8-PSL1 complex and SIX13^22-293 with and without Kex2 protease using 150 nL protein solution and 150 nL reservoir solution sitting-drop plates at 18°C with commercially available sparse matrix screens.

No crystals were obtained for Avr1^18-242, SIX6^17-225, SIX8^50-141, PSL1^18-111, SIX8-PSL1 complex and SIX13^22-293. Final crystallisation conditions were optimised for Avr3^22-284 (0.2 M ammonium sulfate, 0.1 M Bis-Tris pH 6.5, 25% (w/v) PEG 3350), Avr1^59-242 (0.2 M ammonium sulfate, 0.1 M sodium acetate pH 4.5, 17.5% (w/v) PEG 4000), SIX6^58-225 (0.2 M ammonium tartrate and 20% (w/v) PEG 3350), SIX8_C58S^50-141 (0.17 M ammonium sulfate, 15% (v/v) glycerol and 25.5% (w/v) PEG 4000), SIX13 (0.2 M lithium sulfate, 0.1 M Bis-Tris pH 6.5, 25% (w/v) PEG 3350) and PSL1_C37S^18-111 (70% (v/v) MPD and 0.1 M HEPES pH 7.5). Detailed crystallisation optimisation can be found in the supplementary methods.

Crystals were transferred into a cryoprotectant solution containing reservoir solution and 15% (v/v) ethylene glycol, 20% (v/v) glycerol or 10% (v/v) ethylene glycol and 10% (v/v) glycerol. No cryoprotecting was required for SIX8_C58S^50-141 and PSL1_C37S^18-111 crystals as the conditions contained sufficient cryoprotectant (glycerol and MPD, respectively) within the crystallisation condition. For experimental phasing, Avr3^22-284 and Avr1^59-242 crystals were soaked in a cryoprotectant solution containing 0.5 M or 1 M sodium bromide and vitrified in liquid nitrogen. The datasets for bromide-soaked crystals were collected on the MX1 beamline at the Australian Synchrotron [82] (S1 Table). The datasets were processed in XDS [83] and scaled with Aimless in the CCP4 suite [84, 85]. The CRANK2 pipeline in CCP4 was used for bromide-based SAD phasing [86, 87]. Models were refined using phenix.refine in the PHENIX package [88] and model building between refinement rounds was done in COOT [89]. The models were used as a template for molecular replacement against high resolution native datasets collected on the MX2 beamline at the Australian Synchrotron [90]. Automatic model building was done using AutoBuild [91], and subsequent models were refined with phenix.refine and COOT. For SIX6^58-225 and SIX8_C58S^50-141, high confidence ab initio models were generated with AlphaFold2 (S3 Fig), which was used as a template for molecular replacement against a native dataset collected on the MX2 beamline at the Australian Synchrotron. The resultant structure was refined as described above.

Structural modelling and structural alignment

Structural models were generated with Google DeepMind’s AlphaFold2 using the amino acid sequences of SIX effectors and candidates without the signal peptide, as predicted by SignalP-5.0 [92] and predicted pro-domain by searching for a Kex2 cleavage motif (KR, RR or LxxR) if present [30] (S3 Table; S6 Fig). For AlphaFold2 predictions the full databases were used for multiple sequence alignment (MSA) construction. All templates downloaded on July 20, 2021 were allowed for structural modelling. For each of the proteins, we produced five models and selected the best model (ranked_0.pdb). Pairwise alignments of the structural models generated by AlphaFold2 and the experimentally determined structures of Avr1 (PDB code: 7T6A), Avr3 (PDB code: 7T69), SIX6 (PDB code: 8EBB) and SIX8 (PDB code: 8EB9) were generated using the DALI server all against all function [33]. Structural similarity between the pairwise alignments were measured using Z-scores from the DALI server.

Distribution of FOLD family members across fungi

Structure based searches to determine the distribution of FOLD effectors across other phytopathogens was carried out by searching the experimentally determined Avr1, Avr3 and SIX6 structures against available structure databases (Uniprot50, Proteome, Swiss-Prot) using the Foldseek webserver [41] using a 3Di search limited to fungi. An e-value cut off of 0.01 was used, and non-plant associated fungi were removed as well as duplicated results for final analysis. Proteins below 100 amino acids, and above 500 amino acids were filtered out and remaining structural hits were manually inspected for similarity to FOLD effectors.

Interaction studies between PSL1 and SIX8

To investigate the PSL1 and SIX8 interaction in vitro, ∼140 µg of PSL1^18-111 and SIX8^50-141 individually, and ∼140 µg PSL1^18-111 and 140 µg of SIX8^50-141 together were injected onto a Superdex 75 Increase 10/300 (Cytiva) column pre-equilibrated in 20 mM HEPES pH 7.5, 150 mM NaCl, after a 30 min room temperature incubation. To investigate the residues responsible for the interaction, SIX8_C58S^50-141 and PSL1_C37S^18-111 mutants were used instead. Samples across the peaks were then analysed by Coomassie-stained SDS-PAGE. To investigate the mode of interaction, PSL1 and SIX8 proteins and mutants at 10 µM were incubated individually or together for 1 hour at room temperature. An unrelated protein with a free cysteine (AvrSr50^RKQQC) [59] was used to assess the specificity of the PSL1-SIX8 interaction. Proteins were analysed by intact mass spectrometry with or without the addition of DTT as described above.

Agrobacterium-mediated gene expression in N. benthamiana

Agrobacterium tumefaciens cultures containing the pSL constructs and pSOUP [93] were diluted to an OD600 of 1.0 in 10 mM MES pH 5.5 buffer containing 10 mM MgCl2 and 200 μM acetosyringone and incubated in the dark for 2 hours. For co-infiltrations, cultures were mixed together in equal volumes. Resuspensions were infiltrated into 4 – 5-week-old N. benthamiana leaves. Infiltrated plants were kept in a 25°C controlled temperature room with a 16-hour photoperiod. Leaves were imaged 4 – 7 dpi.

Ion leakage assay

Six biological replicates each consisting of three leaf discs (7 mm diameter) were harvested from leaves infiltrated with Agrobacterium 20 - 24 hours post infiltration and incubated in 7 mL of water in a 6-well culture plate. The water was replaced after 40 - 60 min. The leaf discs were incubated in water at room temperature and the conductivity was measured after 24 - 48 hours.

Immunoblot analysis of proteins expressed in N. benthamiana

N. benthamiana leaves infiltrated with A. tumefaciens cultures were harvested 3 dpi. Leaf tissue was frozen in liquid nitrogen, ground into a powder and resuspended in 3x Laemmli buffer containing 0.2 mM DTT and 5 M urea to extract proteins. Samples were boiled for 10 min and centrifuged at 13000 xg to remove leaf debris. Proteins were separated by SDS-PAGE and transferred by electroblotting onto PVDF membranes. Protein blots were probed with anti-HA antibodies conjugated to horseradish peroxidase (Roche, Switzerland, 12013819001, 1:4000). Immunoblots were visualised with Pierce ECL Plus Western Blotting Substrate (Thermo Fisher Scientific) as described by the manufacturer. Membranes were stained with Ponceau S to assess protein loading.

Acknowledgements

This work was supported by the Australian Research Council (ARC DP200100388 D.J./S.W.) and the Australian Academy of Science (Thomas Davies Grant). S.W. was funded by an ARC Future Fellowship (FT200100135) and supported by the ANU Future Scheme (35665). L.M. was funded by an ARC Discovery Early Career Researcher Award (DE170101165). A.S. was a recipient of the AINSE Honours Scholarship Program, and D.Y. and C.M. held an AINSE Postgraduate Research Award. P.K. was supported by a Netaji Subhas ICAR International Fellowship. The authors acknowledge the use of the ANU crystallisation facility. This research was undertaken in part using the MX2 beamline at the Australian Synchrotron, part of ANSTO, and made use of the Australian Cancer Research Foundation (ACRF) detector. The authors acknowledge use of the Australian Synchrotron MX facility and thank the staff for their support. The coordinates and structure factors for Avr1, Avr3, SIX6 and SIX8 have been deposited in the Protein Data Bank with accession number 7T6A, 7T69, 8EBB and 8EB9, respectively.

Supplementary Figures

S1 Fig.

S2 Fig.

S3 Fig.

S4 Fig.

S5 Fig.

S6 Fig.

S7 Fig.

S8 Fig.

S9 Fig.

S10 Fig.

S11 Fig.

S12 Fig.

S13 Fig.

S14 Fig.

Supplementary tables

S1 Table.

S3 Table.

S4 Table.

S5 Table.