Sparse generative modeling via parameter reduction of Boltzmann machines: Application to protein-sequence families

Boltzmann machines (BMs) are widely used as generative models. For example, pairwise Potts models (PMs), which are instances of the BM class, provide accurate statistical models of families of evolutionarily related protein sequences. Their parameters are the local fields, which describe site-specific patterns of amino acid conservation, and the two-site couplings, which mirror the coevolution between pairs of sites. This coevolution reflects structural and functional constraints acting on protein sequences during evolution. The most conservative choice to describe the coevolution signal is to include all possible two-site couplings into the PM. This choice, typical of what is known as Direct Coupling Analysis, has been successful for predicting residue contacts in the three-dimensional structure, mutational effects, and generating new functional sequences. However, the resulting PM suffers from important overfitting effects: many couplings are small, noisy, and hardly interpretable; the PM is close to a critical point, meaning that it is highly sensitive to small parameter perturbations. In this work, we introduce a general parameter-reduction procedure for BMs, via a controlled iterative decimation of the less statistically significant couplings, identified by an information-based criterion that selects either weak or statistically unsupported couplings. For several protein families, our procedure allows one to remove more than $90\%$ of the PM couplings, while preserving the predictive and generative properties of the original dense PM, and the resulting model is far away from criticality, hence more robust to noise.

Figure 1

Fitting and generative quality for PF00076: Pearson correlation coefficient between model and data frequencies as a function of the model density. The one-site frequencies $f_i\left(a\right)$ are directly fitted. The two-site connected correlations $c_{ij}(a,b)$ are fully fitted by the densest model, while only a fraction of them are fitted for the sparse models at $d<1$. The three-site connected correlations $c_{ijk}(a,b,c)$ are never fitted. The generative performance of the model is essentially unchanged down to a density of 10% and slowly decays for even sparser couplings. However, even down to $d=1\%$, the Pearson coefficients remain at remarkably high values above 95% for the two-site correlations, and above 84% for the three-site correlations.

Figure 2

Contact prediction for PF00076: Positive predictive values (PPVs) for several model densities, i.e., the fraction of true positives among the highest-ranking $k$ pairs $(i,j)$ of sites, when ordered by decreasing $F_{ij}^{\mathrm{APC}}$. Even the most sparse model, with only 1.6% of couplings, shows an excellent performance at contact prediction. The curve for plmDCA, a standard DCA approach for contact prediction, is shown for reference and gives comparable results.

Figure 3

Coupling distributions: Distribution of couplings corresponding to true contacts (top) and to noncontacts (bottom) in the three-dimensional protein fold, for the initial PM with density $d_{\max }=91\%$ and a sparser model having density $d=7\%$ associated with a reasonably accurate contact prediction.

Figure 4

Criticality: Heat capacity as a function of temperature for models with different density. The densest models show a strong peak of specific heat close to the reference scale $T=1$, which is a signature of criticality: the model is extremely sensitive to a small change of couplings, due to overfitting. On the other hand, sparse models display a much smaller peak, which is also shifted away from $T=1$ towards lower temperatures, indicating a better robustness of the learning.

Figure 5

Single and double mutations: Spearman correlation between the experimental fitnesses and the model predictions as a function of the model density, for both single and double mutants of the PABP, a member of the PF00076 family. The dashed lines show the same correlations for a profile model ($d=0$) as a reference.

Figure 6

Reconstruction performances measured by the true positive rate (a), false negative rate (b), false positive rate (c), true negative rate (d), and ${\ell }_2$-norm between the true and inferred parameters (f), as a function of the density for the controlled experiment over a Viana-Bray model. We use different colors according to the value of $M$, the number of sequences used within the learning step of the fully connected and of the sparse model. The vertical line represents the density of the true model $d_{\mathrm{true}}=0.0816$. In (e) we show the scatter plot of the true couplings against the parameters learned by the $M=10000$ run for $d=d_{\mathrm{true}}$.

Figure 7

Averaged Hamming distances between an equilibrium sequence at time $t=0$ and the evolved sequence after $t$ MC sweeps for (a) PF00076 (b) PF00014 (c) PF00072 (d) PF00595, and (e) PF13354. The average is computed using ${10}^4$ independent MC chains.

Figure 8

Pearson correlation coefficients between the three chosen metrics (first moments, two-site, and three-site connected correlations) of the data and the sparse models as a function of the density. Each panel (a), (b), (c), and (d) is associated with a different family, respectively, PF00014, PF00072, PF00595, and PF13354.

Figure 9

Positive predictive value (PPV) curve for (a) PF00014, (b) PF00072, (c) PF00595, and (d) PF13354 associated with the contact prediction of several sparse models, from yellow to black lines. As a comparison we show the PPV curve (red line) obtained by the state-of-the-art method for this task, plmDCA.

Figure 10

Distribution of the couplings associated with residues physically in contact (labeled “1” histograms) and with residues not in contact (labeled as “2” histograms) for two different densities. Panels (a), (b), (c), and (d) refer to PF00014, PF00072, PF00595, and PF13354, respectively.

Figure 11

Heat capacity as a function of the temperature $T$ for the other protein families PF00014, PF00072, PF00595, and PF13354 in panels (a), (b), (c), and (d) respectively.

Figure 12

Spearman correlation coefficient between the energy variations computed according to the sparse models (spBM lines) and the experimentally determined fitness variations of a set of single mutants, for the TEM1 (PF13354) and PSD95 (PF00595) proteins. The dashed lines show the results of the Spearman correlation coefficients when the energy variations are computed with the profile models of the corresponding families.

Figure 13

Left panels: Pearson correlation coefficients between (a) the first moments, and (b) the two- and (c) three-site connected correlations of each model (varying the density) compared to the empirical data. The lines are colored according to the metrics used within the decimation procedure: the symmetric Kullback-Leibler distance, and the strength of the couplings or the two-site frequencies. Right panels: Pearson correlation coefficients between (d) the first moments and (e) the two- and (f) three-site connected correlations varying the initial condition of the dense model learning. All data are for the PF00076 family.

Figure 14

Positive predictive value for each decimation procedure (a) and for each initial condition (b) for the PF00076 family.

Figure 15

Pearson correlations (a) and positive predictive value (b) for the decimation via ${\ell }_1$-norm regularization, for the PF00076 family, compared with those reported in the main text.

Figure 16

Cumulative density function of the logarithms of the two-site frequencies associated with the decimated couplings for the sparse models having densities $d=90.7\%, d=69.9\%, d=49.7\%, d=12.2\%$, and $d=3.2\%$. The data refer to the PF00076, and the decimation is performed according to the standard protocol described in the main text.

Figure 17

Pearson correlation coefficients for the online learning (blue line), for $K=10,20,40$, compared to the converged run (red line) as a function of the density, for the one-site frequencies (a), two-site (b), and three-site (c) connected correlations.

Figure 18

Comparison between the PPV curve obtained for the online (for $K=10,20,40$) and converged runs at two different densities, 72% (a) and 3.2% (b). All data are for the PF00076 family.

Figure 19

Sequence variability for all the protein families considered in this work. We plot $D_{\mathrm{ts}}, D_{\mathrm{st}}$, and $D_{\mathrm{ss}}$ as a function of the density using blue, orange, and green lines, respectively, for PF00076 (a), PF00014 (b), PF00072 (c), PF00595 (d), and PF13354 (e).