A large accessory protein interactome is rewired across environments

Abstract

To characterize how protein-protein interaction (PPI) networks change, we quantified the relative PPI abundance of 1.6 million protein pairs in the yeast Saccharomyces cerevisiae across nine growth conditions, with replication, for a total of 44 million measurements. Our multi-condition screen identified 13,764 pairwise PPIs, a threefold increase over PPIs identified in one condition. A few 'immutable' PPIs are present across all conditions, while most 'mutable' PPIs are rarely observed. Immutable PPIs aggregate into highly connected 'core' network modules, with most network remodeling occurring within a loosely connected 'accessory' module. Mutable PPIs are less likely to co-express, co-localize, and be explained by simple mass action kinetics, and more likely to contain proteins with intrinsically disordered regions, implying that environment-dependent association and binding is critical to cellular adaptation. Our results show that protein interactomes are larger than previously thought and contain highly dynamic regions that reorganize to drive or respond to cellular changes.

Keywords: S. cerevisiae; computational biology; dynamic ppis; immutable PPIs; mutable PPIs; network rewiring; protein interactome; protein-protein interactions; systems biology.

Figure 1.. PPiSeq.

(A) A cartoon of PPiSeq yeast library construction. Strains from the protein interactome collection are individually mated to strains from the double barcoder collection and sporulated to recover haploids that contain a mDHFR-tagged protein and a barcode. Haploids are mated as pools. In diploids, expression of Cre recombinase causes recombination between homologous chromosomes at the loxP locus, resulting in a contiguous double barcode that marks the mDHFR-tagged protein pair. (B) Representative double barcode frequency trajectories over twelve generations of competitive growth. Trajectories are used to calculate a quantitative fitness for each double barcoded strain. (C) Standard error of fitness estimates of protein pairs. The blue and red lines represent the median standard error for a sliding window (width = 0.05) of all fitness-ranked protein pairs and of only the positive protein-protein interactions, respectively. (D) Estimated fitness of strains with different double barcodes representing the same protein pair in the same pooled growth. Positive protein pairs are randomly selected within a fitness window. ORF x Null is a violin plot of the fitness distribution of all interactions with a mDHFR fragment that is not tethered to a yeast protein. DHFR(-) is yeast strains that lack any mDHFR fragment. DHFR(+) is yeast strains that contain a full length mDHFR under a strong promoter. (E) Density plot of the fitness of double barcodes that represent the same putative PPI in the same pooled growth. In (B–E), the data in SD environment are used. (F) Density plot of the normalized mean fitness of the same PPI between two pooled growth cultures in SD environment. PPIs detected in either one growth culture are included. (G) Venn diagram of the number of PPIs identified within our search space by PPiSeq in nine environments (magenta), PPiSeq in SD environment (pink), the interactome-scale protein-fragment complementation screen (PCA, yellow), and the BioGRID database excluding any PPIs previously detected by PCA (blue).

Figure 1—figure supplement 1.. Double barcodes and protein pairs in the PPiSeq library.

(A) Distribution of the initial double barcode count of the PPiSeq library in SD environment at a sequencing depth of 209,899,687 reads. (B) Number of barcodes per protein pair in the PPiSeq library. Spike-in control protein pairs are not included in the plot.

Figure 1—figure supplement 2.. Standard error of fitness estimates of protein pairs in each environment.

The blue and red lines represent the median standard error for a sliding window (width = 0.05) of all fitness ranked protein pairs and of only the positive protein-protein interactions, respectively.

Figure 1—figure supplement 3.. Density plot of the fitness of double barcodes that represent the same positive PPI in the same pooled growth of each environment.

Figure 1—figure supplement 4.. Comparison of PPiSeq data in SD condition to other PPI datasets.

(A) Heatmap of overlap coefficients across different datasets. PPiSeq SD1 and PPiSeq SD2 are PPIs identified from two replicate growth cultures. PPiSeq SD-merge is PPIs identfied from the merged data of two replicate growth cultures. PCA is PPIs identified from mDHFR-PCA colony screening (Tarassov et al., 2008). Y2H is PPIs identifed from the latest large-scale yeast-two-hybrid screen (Yu et al., 2008). PRS is a newly constructed ‘bronze standard’ high-confidence positive reference set. Only PPIs within our search space are considered and the numbers are shown under each dataset. (B) Venn diagram of the number of PPIs identified in the two SD PPiSeq replicates and in each of PCA, PRS, or Y2H. (C) Barplot of the overlap coefficients between PPiSeq and other methods for PPIs that do or do not across the two SD PPiSeq replicates.

Figure 1—figure supplement 5.. The OD 600 trajectories of DHFR(-) strain in various conditions with and without 0.5 μg/mL methotrexate.

Figure 2.. Functional enrichment of PPIs detected by PPiSeq.

PPI enrichment (red) and variability (blue) across environments of gene ontology cellular compartments (A) and biological processes (B). Red node size is the percent of interacting protein pairs (interaction density) observed for a given pair of GO terms and the node color is the p-value of this percent over a random expectation. Blue node size and color are the variability (coefficient of variation, CV) of interaction densities across nine environments tested. GO terms are hierarchically clustered by the interaction density (red dots). Boxes mark frequently interacting and invariable cellular compartments and biological processes involved in membrane transport and protein maturation (blue and green) and cell division (purple). Barplots show the mean CV of interaction densities for each GO term across all other GO terms. Orange, black and brown triangles highlight three different groups of related GO terms: chromosome, transcription, and translation, respectively.

Figure 3.. A large accessory protein interactome.

(A) Barplot of PPI number binned by the number of environments in which a PPI is observed. Colors indicate PPIs called by both PPiSeq and BioGRID inclusive of mDHFR-PCA (red), PPIs called by PPiSeq that scored high but were not called by mDHFR-PCA (yellow), and PPIs called by PPiSeq that scored low by mDHFR-PCA (blue). (B) Validation rates of PPIs binned by the number of environments in which a PPI is observed. Validations were performed using OD600 trajectories of clones grown in multi-well plates. (C) Mutable and less mutable PPIs form distinct modules in the network. PPIs that are detected in at least five environments (red edges) form two tight core modules. PPIs that are detected in fewer than five environments (blue edges) form a less connected accessory module. Proteins in different modules are labeled with different shapes and colors. The network uses an edge-weighted spring embedded layout. (D) Number of PPIs within and between each community. PPIs detected in at least or fewer than five environments are shown in red and blue, respectively. The size of the square or circle is proportional to the number of PPIs. The number below each community is the number of proteins within each community. (E) Scatter plot of degrees and mutability scores of proteins in each community.

Figure 3—figure supplement 1.. PPIs across conditions.

(A) Barplot of number of PPIs in each environment binned by the number of environments in which a PPI is observed. (B) Heatmap of fitness values of all detected PPIs across different environments. PPIs (rows) and environments (columns) are hierarchically clustered by the fitness values across environments. (C) Scatter plot of mean fitness values of the same PPI across two different growth conditions. Colors indicate in which condition(s) PPIs are called. PPIs that have been detected in at least one environment are shown in (B) and (C). Negative values and missing measurements are replaced with zeros.

Figure 3—figure supplement 2.. Validating PPIs.

(A) Boxplots and univariate scatterplots of fitness values of PPIs binned by the number of environments in which a PPI is observed. The bottom of each box, the line drawn in the box, and the top of the box represent the 1st, 2nd, and 3rd quartiles, respectively. The whiskers extend to 1.5 times the interquartile range (from the 1st to 3rd quartile). The fitness for each PPI is calculated by taking the mean of the fitness values for all environments where that PPI was detected. (B) Validation rates of PPIs binned by the number of environments in which a PPI is observed. Validations were performed using OD600 trajectories of clones grown in multi-well plates. PPIs that have been previously reported in BioGRID (red) or are previously unreported (blue) are shown. (C) Validation of 20 randomly chosen PPIs that were only detected in SD by PPiSeq. Validations use OD600 trajectories of clones grown in different environments. Red and blue boxes represent positive and negative PPI detection, respectively. (D) Density plot of the increase in the relative area under the growth curve (AUC) against a negative control strain for the 20 PPIs shown in (C). Dashed vertical lines represent the mean AUC increase for an environment.

Figure 3—figure supplement 3.. Predicting validation rates.

(A and B) Distributions of mean fitness value (f) and number of environments in which the PPI is detected (n) for validated and unvalidated PPIs. (C) Comparison between observed validation rates and predicted validation rates for 502 PPIs binned by n. Predicted validation rates were calculated using the mean f and n within each bin. (D) Barplot of the PPI number and the predicted true positive PPI number binned by the number of environments in which a PPI is observed.

Figure 3—figure supplement 4.. Mutable PPIs outnumber immutable PPIs in higher confidence PPI networks.

Barplot of the PPI number binned by the number of environments in which a PPI is observed in the multi-condition network made by either higher confidence PPI calls (A) or by excluding the 16°C condition (B). Colors indicate PPIs called by both PPiSeq and BioGRID inclusive of mDHFR-PCA (red), PPIs called by PPiSeq that scored high but were not called by mDHFR-PCA (yellow), and PPIs called by PPiSeq that scored low and were not called by mDHFR-PCA (blue).

Figure 3—figure supplement 5.. PPIs with a similar mutability are more likely to be connected.

(A) Degree density of proteins binned by the number of environments in which a PPI is detected. (B) Degree density of all proteins that are neighbors of proteins binned by the number of environments in which a PPI is detected. (C) Mutability score density of proteins binned by the number of environments in which a PPI is detected. (D) Mutability score density of all proteins that are neighbors of proteins binned by the number of environments in which a PPI is detected. Any unique protein that participates in a PPI within each bin was counted. The degree (A) and mutability score (C) of a protein were obtained from a multi-environment network that includes all PPIs detected in at least one environment. The degree (B) and mutability score (D) for a target protein’s neighbor was calculated as above, only the interaction between the target protein and its neighbor was first removed.

Figure 3—figure supplement 6.. The multi-environment PPI network contains three major communities with different mutability scores.

Boxplots and univariate scatterplots of mutability scores of communities with at least 10 proteins identified by (A) Fast-Greedy (B) Walktrap, and (C) InfoMAP algorithms. The bottom of each box, the line drawn in the box, and the top of the box represent the 1st, 2nd, and 3rd quartiles, respectively. The whiskers extend to 1.5 times the interquartile range.

Figure 4.. Properties of mutable and less mutable PPIs.

(A) The co-expression mutual rank for PPIs binned by the number of environments in which the PPI is detected. A higher mutual rank means worse co-expression. Notches are the 95% confidence interval for the median, hinges correspond to the first and third quartiles, and whiskers extend 1.5 times the interquartile range. (B) The percent of protein pairs that have been found colocalized by gene ontology (GO Slim, dashed line) and fluorescence (solid line) (Chong et al., 2015). (C) Spearman correlation between the PPI mutability score and other gene features, binned a gene’s PPI degree. In (B) and (C), the error bars are the standard deviation from 1000 bootstrapped data sets. (D) Pearson correlation between a PPI’s fitness and geometric mean abundance of two interacting proteins in Ho et al., 2018, binned by the number of environments in which a PPI is detected. (E) Examples of non-significant (Erv25 x Shr3) and significant (Akr1 x Any1) predictions. Observed heterodimer fitness (S_AB) is plotted against the expectation based on the geometric mean of the two constituent homodimer fitnesses (S_AA and S_BB). (F) Percent of heterodimers whose fitness changes can be significantly predicted by the geometric mean of the two constituent homodimers, binned by the number of environments in which a PPI is observed.

Figure 4—figure supplement 1.. The co-expression mutual rank for PPIs detected in each condition binned by the number of environments in which the PPI is detected.

Notches are the 95% confidence interval for the median, hinges correspond to the first and third quartiles, and whiskers extend 1.5 times the interquartile range.

Figure 4—figure supplement 2.. Mutable PPIs and their properties for higher confidence PPI calls.

(A) The co-expression mutual rank for PPIs binned by the number of environments in which the PPI is detected. A higher mutual rank means worse co-expression. Notches are the 95% confidence interval for the median, hinges correspond to the first and third quartiles, and whiskers extend 1.5 times the interquartile range. (B) The percent of protein pairs that have been found colocalized by gene ontology (GO Slim, dashed line) and fluorescence (solid line) (Chong et al., 2015). (C) Spearman correlation between the protein’s mutability score and other gene features. (D) Spearman correlation between the PPI mutability score and other gene features, binned a gene’s PPI degree. In (B–D), the error bars are the standard deviation from 1000 bootstrapped data sets.

Figure 4—figure supplement 3.. Mutable PPIs and their properties, excluding the 16°C condition.

(A) The co-expression mutual rank for PPIs binned by the number of environments in which the PPI is detected. A higher mutual rank means worse co-expression. Notches are the 95% confidence interval for the median, hinges correspond to the first and third quartiles, and whiskers extend 1.5 times the interquartile range. (B) The percent of protein pairs that have been found colocalized by gene ontology (GO Slim, dashed line) and fluorescence (solid line) (Chong et al., 2015). (C) Spearman correlation between the protein’s mutability score and other gene features. (D) Spearman correlation between the PPI mutability score and other gene features, binned by a gene’s PPI degree. In (B–D), the error bars are the standard deviation from 1000 bootstrapped data sets.

Figure 4—figure supplement 4.. Spearman correlation between the protein’s mutability score and other gene features.

The error bars are the standard deviation from 1000 bootstrapped data sets.

Figure 4—figure supplement 5.. Proteins that participate in multiple complexes are distributed over a wide range of complexes.

(A) Barplot of the number of proteins binned by the number of protein complexes in which a protein participates in. Color represents a protein’s PPI degree. (B) Histogram of number of genes per protein complex (Costanzo et al., 2016). (C) Histogram of number of genes per protein complex for genes included in the PPiSeq screen and involved in at least two protein complexes (Costanzo et al., 2016).

Figure 4—figure supplement 6.. Exploring the relationship of protein abundance, PPI abundance, and PPI mutability.

(A) Density plot of the fitness of a PPI strain in the SD environment against the geometric mean protein abundance from Ho et al., 2018. Colors are the density of points in hexagonal bins. (B) Linear fits to a mass-action kinetic model where the x-axis is the heterodimer PPI fitness expected from the homodimer fitnesses of the constituent proteins, and the y-axis is the measured heterodimer fitness. The left panel contains heterodimers significantly explained by the model that do not require a significant intercept (FDR < 0.05, see Materials and methods). Colors are the R² of the mass-action kinetics model fit. (C) The coefficients of each scaled feature in a logistic model predicting a good fit to the mass-action kinetics model, as fit by ‘glm’ function in R. (D) Density plot of the geometric mean abundance of a heterodimer pair against the Pearson correlation between the predicted and observed heterodimer fitness across conditions. Colors are the density of points in hexagonal bins. Red line is a Deming regression, r is the Pearson correlation. (E) Explained PPIs are composed of less abundant proteins. Box and dot plot of the mean protein abundance of a heterodimer for PPIs that are explained and not explained by the mass-action kinetics model. Boxplot summarizes the first, second, and third quartiles. ***p<10⁻⁹ Wilcoxon signed-rank test.

Figure 5.. Carbohydrate transport network rewiring as captured by PPiSeq.

(A) Heatmap of abundances (fitnesses) of PPIs involved in carbohydrate transport across different environments. (B) Boxplots of fitnesses of PPIs involving Hxt proteins in SD, Raffinose and NaCl environments. The bottom of each box, the line drawn in the box, and the top of the box represent the 1st, 2nd, and 3rd quartiles, respectively. The whiskers extend to ±1.5 times the interquartile range. (C) Circular network plots of PPIs containing Hxt proteins in SD, Raffinose, and NaCl environments. Nodes are proteins and colors are as in (A). Node size is proportional to its degree in the multi-environment PPI network. Edge width is proportional to abundance in each environment. (D) Scatter plot of fitness changes relative to SD as measured by PPiSeq and clonal growth dynamics for randomly chosen carbohydrate-transport PPIs in Raffinose (80 PPIs) and NaCl (90 PPIs).

Figure 6.. The estimated number of true PPIs discovered by PPiSeq using repeated sampling of data in permuted orders of environment addition.

Boxplots summarize the distribution of the number of unique PPIs across permutations. The bottom of each box, the line drawn in the box, and the top of the box represent the 1st, 2nd, and 3rd quartiles, respectively. The whiskers extend to ±1.5 times the interquartile range. Overlayed solid red lines and dashed red lines are the Kindt exact accumulation curves and the bootstrap estimators of the total number of unique PPIs across infinite environments for each simulation, respectively.

Appendix 1—figure 1.. Defining a dynamic threshold for PPI calling.

(A) A discrete combination of a fitness threshold (f) and a p-value threshold (p) results in a PPV. Colored lines are fitness and p-value thresholds that result in the same PPV in SD. (B) Density plot of all f and p combinations that result in a PPV of 0.7 using 50 different random reference sets in SD. The black line is the fitted sigmoid model that is used for the dynamic threshold.

Appendix 1—figure 2.. Density plots of the dynamic thresholds in each environment.

Data were split into two groups: p < −4 and p >= −4. For p >= −4, as in Appendix 1—figure 1B, a sigmoidal function was fit to f and p combinations that result in the same PPV value. For p < −4, the fitness threshold was set to equal the minimum fitness value when p >= −4. The dynamic threshold that results in the maximum MCC in each environment was shown in the plot.

Appendix 1—figure 3.. The precision-recall curves of dynamic thresholds in each environment.

Red asterisks mark the thresholds with maximal Matthews correlation coefficients.

Appendix 1—figure 4.. Dynamic thresholds (red) of f and p have a higher positive predictive value (PPV) than most discrete combinations (blue).

Points represent the PPVs for dynamic thresholding and for all combinations of discrete fitness and p-value thresholds underlying a constant range of false positive rates obtained from the optimal dynamic threshold in each environment.