MultiOmics integration
Following the holo-omics approach, we integrated the metabolomics and root metabarcoding datasets using the DIABLO framework (Singh et al. 2019) with multiblock sparse PLS-DA (partial least squares discriminant analysis). Model tuning helped us to select features from metabarcoding and metabolomics to improve the modelling of differences between groups. Noteworthy, the datasets are highly correlated for all three components (Figure 6, components 1 and 2, supplementary figure S5, components 2 and 3). The agreement between metabarcoding and metabolomics is high for all samples and treatments. The combined drought and heat treatment separates from the other three treatments on the first axis, and heat treatment on the second axis, while the drought and control treatments are not well separated by components 1 and 2, but are separated by component 3 (Supplementary figure S5.
As a follow-up analysis to DIABLO, Pearson correlation network analysis was performed. We aimed to identify positive correlations between features assigned to the same treatment by DIABLO. These pairs likely resolve the underlying interaction between the plant metabolome and the rhizosphere microbiota.
Only bacterial genera and plant metabolites significantly assigned to a feature by DIABLO were included in network construction. Three clustering methods (Optimal, Louvain, fast greedy clustering) reproduced four similar modules. The largest was mainly constituted by correlations between features of the control treatment and several heat or drought metabolites. The latter were likely false positives of the DIABLO analyses.
Besides the control module, only metabolites and genera assigned as features of the H+D treatment by DIABLO co-occurred in a single module. Thus, crosslinks between metabolite and bacterial datasets were characteristic of the H+D treatment appearing as chords within Figure 7 . Additionally, a fraction of drought-associated bacterial genera and heat-associated metabolites formed connections within this majorly H+D-associated module. The module consisted of 24 nodes in total (Table 1). Among them were predominantly Nitrogen-Containing Secondary Compounds for biosynthesis, Phenylpropanoid Derivative biosynthesis, Polyketide biosynthesis and Proteobacteria .
In contrast, the majority of drought-associated metabolites or heat-associated bacterial genera formed no crosslinks between the same feature category across the metabolite and metabarcoding datasets. Instead, drought-associated metabolites correlated with other metabolites and heat-associated bacterial genera interacted with other bacteria, as indicated by characteristic arcs within Figure 7.
In conclusion, a small core of Proteobacteria was associated with Nitrogen-Containing Secondary Compounds in response to combined H+D stress. The combined treatment appears to further select plant metabolism responsive to heat metabolism, but rhizosphere bacteria responsive to D.