Estimating statistical significance of local protein profile-profile alignments

Background Alignment of sequence families described by profiles provides a sensitive means for establishing homology between proteins and is important in protein evolutionary, structural, and functional studies. In the context of a steadily growing amount of sequence data, estimating the statistical significance of alignments, including profile-profile alignments, plays a key role in alignment-based homology search algorithms. Still, it is an open question as to what and whether one type of distribution governs profile-profile alignment score, especially when profile-profile substitution scores involve such terms as secondary structure predictions. Results This study presents a methodology for estimating the statistical significance of this type of alignments. The methodology rests on a new algorithm developed for generating random profiles such that their alignment scores are distributed similarly to those obtained for real unrelated profiles. We show that improvements in statistical accuracy and sensitivity and high-quality alignment rate result from statistically characterizing alignments by establishing the dependence of statistical parameters on various measures associated with both individual and pairwise profile characteristics. Implemented in the COMER software, the proposed methodology yielded an increase of up to 34.2% in the number of true positives and up to 61.8% in the number of high-quality alignments with respect to the previous version of the COMER method. Conclusions The more accurate estimation of statistical significance is implemented in the COMER method, which is now more sensitive and provides an increased rate of high-quality profile-profile alignments. The results of the present study also suggest directions for future research. Electronic supplementary material The online version of this article (10.1186/s12859-019-2913-3) contains supplementary material, which is available to authorized users.


Introduction
Establishing homology between proteins is essential in evolutionary and highly important in structural and functional studies of proteins and their complexes. Alignment, in general, represents the most fundamental way to deduce homology directly or indirectly. Sequence alignment has proved indispensable in annotating uncharacterized proteins and paved the way for alignment of sequence families described by profiles, constituting the basis for inferring protein structure and function by homology.
In the presence of a large and steadily growing amount of sequence data, estimating the statistical significance of alignments plays a prominent role in alignment-based homology search algorithms (Karlin, 2005). For an alignment with a particular similarity score, it provides a probability or related quantity indicating how likely a chance alignment with the same or greater score is to be observed.
The limiting distribution of the ungapped local alignment score S (nonlattice) for large sequence lengths m * mindaugas.margelevicius@bti.vu.lt and n has been proved Karlin and Altschul, 1990;Karlin and Brendel, 1992;Dembo et al., 1994) to be the type 1 (Gumbel-type) extreme value distribution (EVD) (Kotz and Nadarajah, 2000) where λ and K (K = e λµ /[mn]; µ, the location parameter under the standard parametrization) are parameters calculable from the score matrix used to align sequences and satisfying the conditions of the negative expected score and existence of at least one positive score. Ample empirical evidence (e.g., Mott, 1992;Altschul and Gish, 1996;Pearson, 1998;Waterman and Vingron, 1994b,a) has suggested that the statistical theory for ungapped alignments applies to gapped alignments, provided gap penalties, or costs, lead to alignment scores in the (local) region of logarithmic growth (Arratia and Waterman, 1994).
This study aims at characterizing the distribution of profile-profile alignment scores to improve statistical accuracy and remote homology detection. The profile achieves its highest potential in recognizing related proteins when its record at each position combines information components obtained by independent numerical methods or sources of information, e.g., secondary structure (SS) assignment, structure/sequence context predictions, contact propensities, and others. We consider this type of profiles, i.e., supplemented with composite data. To the best of our knowledge, studies dedicated to the distribution of scores produced by aligning such profiles have not been reported in the literature, although the characterization of profile-profile alignments assisted at least with SS predictions has attracted much interest.
Our focus is local alignment scores and the specification of their distribution. Hence, we assume that the expected profile-profile alignment score per aligned pair is negative and the probability for a positive score is positive. These assumptions are satisfied in part when the score for a pair of positions of two profiles is (implicitly) log-odds (Karlin, 2005;Meng et al., 2011). Importantly, the score preserves the form of log-odds when it linearly combines different components (e.g., the similarity of two amino acid frequency vectors with that of SSs) as long as each component itself represents a logodds score. Such a construction of composite scores is typical, and, given appropriate gap penalties, we can expect that the assumptions will hold for most profileprofile alignment methods.
Here, we use the COMER method (Margelevičius, 2016) as a means for producing profile-profile alignment scores the algorithms developed in this work employ. The COMER profile at each position encapsulates the transition probabilities of moving to and from the states of insertion and deletion, which for a pair of profiles transform to position-specific gap penalties. Comparing two COMER profiles also includes scoring predicted secondary structure and sequence contexts (Margelevičius, 2018). These properties make the characterization of the distribution of alignment scores a challenging task.
We aim at proposing a general procedure for estimating the statistical significance of profile-profile alignments. Therefore, we do not consider any modifications to the profile model, such as introducing additional HMM states and their transitions for more reliable alignment statistics (Eddy, 2008), which, however, affect alignments and their scores not necessarily leading to increased sensitivity. This fits the procedure's objective to specify an alignment score distribution for profiles representing as much as possible real unrelated proteins.
The importance of the null model of random protein sequences has long been recognized (Mott, 2000;Yu and Hwa, 2001;Metzler, 2006;Eddy, 2008;Meng et al., 2011). Real sequences do not conform to the canonical (and simplest) model, in which the amino acids in the sequence are iid, and exhibit a more complex structure with short-and long-range amino acid correlations. It has been proved that in the limit of infinitely long sequences, ungapped alignment scores for random Markov-dependent sequences are still distributed according to the EVD (Karlin and Dembo, 1992). Gapped alignments of random correlated sequences accounted for by a null model have been shown numerically to correspond to an EVD too (Messer et al., 2007). Notably, the values of the statistical parameters differed substantially from those obtained for iid sequences.
In this work, the null model of protein sequence family does not explicitly incorporate correlations between families and leaves the profile-profile scoring function unchanged. Instead, we randomly generate profiles with properties of real unrelated profiles and take into account the correlations between them by establishing the dependence of statistical parameters on the compositional similarity between the profiles. The null model affects the probability of a score of a pair of profile positions, yet the change in composition does not alter the type of alignment score distribution but values of statistical parameters. The approach proposed here predicts differences in these values.
Overall, this work purposes a procedure for estimating the statistical significance of profile-profile alignments, regarding the effects and factors that influence alignment score distribution, such as edge effects, compositional similarity, and profile attributes (length and the effective number of observations). The procedure builds on an algorithm proposed for generating profiles that resemble real sequence families. The results of the implemented procedure give rise to an analysis of its impact on sensitivity and alignment quality, which we also provide in this paper.

Methods
The concept of estimating the statistical significance of profile-profile alignments, as proposed here, is illustrated in Figure 1. The significance of the alignment between a profile of length l1 and effective number of observations (ENO) n1 and another profile of length l2 and ENO n2 follows from combining the significance of the alignment score and that of the number of positive substitution scores in the alignment. The distribution of alignment scores depends on both the length of profiles being compared and their ENOs (Sadreyev and Grishin, 2008) (Supplementary Section S2.1 defines ENO). It also depends on the compositional similarity of profiles. Therefore, the significance of the alignment score of profiles 1 and 2 is calculated using statistics obtained from aligning unrelated profiles of the same length and ENO and having the same mutual compositional similarity.
The questions related to the proposed concept are the following: 1. How is the compositional similarity between profiles calculated? 2. What types of distributions govern alignment score and the number of positive substitution scores? 3. How are unrelated profiles generated? 4. How are statistics that are obtained from different simulations combined? 5. How are statistics derived for any values of length, ENO, and compositional similarity? 6. How are dependent p-values combined?
We consider them in detail below. More details can be found in Supplementary Section S4.

Compositional similarity
The statistical parameter λu ≡ λ of the limiting distribution of sequence alignment scores (1) has several related meanings (Karlin and Altschul, 1990;Altschul et al., 1997;Schäffer et al., 2001). λu can also be regarded as a measure of compositional similarity (Yu et al., 2006). The value of λu found as the positive solution to X k p(s k ) exp(λus k ) = 1, where {s k } k represent different values of scores in the substitution matrix and p(s k ) is the probability of s k , will decrease as the number of positive substitution scores increases (p(s k ) increases for all k : s k > 0). Therefore, low values of λu may indicate an increased probability for a high-scoring alignment to occur by chance due to compositionally biased regions in two sequences. The parameter λu calculated for a pair of profiles (Sadreyev and Grishin, 2003;Margelevičius and Venclovas, 2010) has the same dependence on composition. We take into account compositional similarity between profiles by specifying the dependence of the distribution of profile-profile alignment scores on it.

Alignment scores of real unrelated profiles
We analyze the distribution of alignment scores of real unrelated profiles (constructed from multiple sequence alignments, MSAs, of real sequences) for two reasons. One is to determine the type of distribution that describes them. The other reason is that these data provide a reference point for generating random profiles whose alignment scores would follow the same type of distribution and be similarly distributed. We found that alignment scores of real unrelated profiles follow an EVD and 85% of goodness-of-fit tests do not reject the null hypothesis. Supplementary Sections S2.3 and S3.1, respectively, describe how real unrelated profiles were obtained and detail the comparison of profiles of different length, ENO, and compositional similarity.

Profile simulation
It is impractical to use real unrelated profiles to produce sufficient data for any values of profile length, ENO, and compositional similarity. Instead, random profiles are generated.
It has been shown that using a null model to generate realistic random profiles improves remote homology detection (Sadreyev and Grishin, 2008). Therefore, the aim is to develop an algorithm for generating random profiles that exhibit properties of real profiles, and the degree of correlation between the random profiles can be controlled.
We achieve that by using a modest number of "seed" profiles constructed for a diverse set of real sequences. Random profiles then result from adding noise to profiles/MSAs generated using these seed profiles as a model. Seeds provide properties of real profiles, whereas noise and randomization represent a means of controlling correlation between random profiles. However, Supplementary Section S3.2 demonstrates the complications that arise when the comparison of generated profiles includes the scoring of SSs predicted for profiles.
Hence, an important feature of the algorithm that we propose (Algorithm 1) is that random profiles result from concatenating fixed-length fragments sampled randomly from a set of generated MSAs regardless of SS predictions the fragments entail. Supplementary Section S4.1 provides additional details, and Supplementary Algorithm S2 defines how MSAs with added noise are produced using a profile model in step 2 of Algorithm 1.
Algorithm 1 Generating random profiles of length l and ENO n given noise level r and fragment length s Input: S sequences. Output: R random profiles.
1. Make profiles (seeds) for S diverse sequences chosen at random for which profiles are obtained to be of a sufficiently large ENO (e.g., 12).
2. Using each of the S seed profiles as a model, generate M MSAs of ENO n with noise r.
3. Using each set of S×M generated MSAs of ENO n as source MSAs, generate a set of R random MSAs of length l in the following way: (a) choose a number j at random uniformly between 1 and S×M ; (b) randomly select a fragment of length s of source MSA j; (c) copy and add the selected fragment to a random MSA being generated; (d) repeat steps (3a)-(3c) until the length of the random MSA becomes l.
4. Construct profiles from the R generated random MSAs.
Algorithm 1 shows that the fragment length s and the level r of noise added to S×M generated MSAs are the parameters that determine the degree of similarity among resulting random profiles. Too large s and/or too small r lead to highly correlated random profiles, while too small s or too large r result in divergent profiles and useless alignment statistics.

Optimizing s and r
Based on the results reported in Section 2.2, we find the optimal s and r using two criteria: (1) the goodness of fit of the EVD to the empirical distribution of alignment scores and (2) the distance between the distribution function obtained for real unrelated profiles and the distribution function obtained from simulations. In this way, profiles generated using the optimal s and r possess features characteristic to real profiles and correlations between them do not dominate.
We calculate the supremum class upper tail Anderson-Darling statistic ADup (Chernobai et al., 2015) (Supplementary Section S5.1) for testing the goodness of fit of the EVD to the data (the first criterion). We normalize it, AD up = ADup/ √ N , by the square root of the number of alignment scores, √ N , to make it independent of sample size. The distance between two empirical distribution functions (the second criterion) is measured by the two-sample Kolmogorov-Smirnov statistic D.
The whole procedure for finding the optimal s and r can be described as follows. 1. Choose values for s and r. 2. Using Algorithm 1, generate random profiles of all ENOs and lengths considered. 3. Align simulated profiles. 4. Fit the EVD in the right tail of alignment score distributions. 5. Calculate the AD up statistic for each alignment score distribution. 6. Train models for predicting the statistical parameters for profiles with given attributes (ENO and length) and compositional similarity. 7. Predict the statistical parameters and calculate the p-value of each alignment of real unrelated profiles. 8. For each combination of pair values of profile ENO and length, calculate the distance D between the empirical distribution function obtained for real unrelated profiles and the distribution function obtained in step 7. 9. Repeat steps 1-8 until the optimal s and r with respect to AD up and D have been found.
We discuss prediction of statistical parameters (step 7) below. Here we note that predictions are used to incorporate compositional similarity into the statistical model of profile-profile alignments.
The results (Figure 2A) obtained using S = 1012 seed profiles (M = 1) suggest that considering a small number of different values suffices to find the optimal s and r: Large s increases correlation between profiles, whereas large r makes them unalignable. Figure 2A shows that the best balance between the two criteria is achieved for s = 9 and r = 0.03.
Sections S6.1 and S6.2 present additional results from simulation experiments. They also show that seeding from three different profiles representing different SCOPe (Fox et al., 2013) classes leads to similar results ( Figure 2B). Using one seed facilitates profile simulation.

Combining statistics
Consider two profiles, one of length l1 and ENO n1 and another one of length l2 and ENO n2. The profiles share compositional similarly λu. Then, the significance of their alignment score is determined from combining statistics obtained (1) from aligning simulated profiles of length l1 and ENO n1 against profiles of length l2 and ENO n2 and (2) from aligning profiles of length l1 against profiles of length l2 with mutual compositional similarly λu, as shown in Figure 1. Here we define this statistical combination. Let A and B index two distributions. The first is obtained for profiles described by the set of attributes {n1, l1; n2, l2} (i.e., profiles of length l1 and ENO n1 have been aligned against profiles of length l2 and ENO n2) and the other for profiles characterized by the set of attributes {λu; l1; l2}. Assume that distributions A and B belong to the family of EVDs. Letμ A andσ A respectively denote the estimates of the location and scale parameters of the EVD corresponding to distribution A. Similarly, letμ B andσ B be the estimates of the statistical parameters of distribution B. Then, the conditional mean estimator of the location parameter given two sets of parameters specifying its distribution in settings A and B iŝ and the corresponding conditional mean estimator of the scale parameter iŝ where a and b are parameters that depend on the parameters specifying the distributions of the location and scale parameters, respectively, in settings A and B (see Supplementary Appendix A for details and a proof).
Statistical significance for profiles of any length l, ENO n, and λu has to be estimated based on the statistics obtained from the observed distributions. As the statistics depend non-linearly on l, n, and λu (Supplementary Section S6), we use low-complexity artificial neural network (NN) models to predict statistical parameters. The models are trained on the estimates obtained (1) from aligning simulated profiles of length l1 ∈ L and ENO n1 ∈ N against profiles of length l2 ∈ L and ENO n2 ∈ N (n2 < n1) and (2) from aligning profiles of length l1 ∈ L against profiles of length l2 ∈ L with λu discretized. We denote the predictions of the location parameter in these two settings byμ A (n1, l1; n2, l2) andμ B (λu; l1; l2), respectively. The notation for the scale parameter is similar. (See also Supplementary Section S4.4.)

Adjustment of predictions
Statistical accuracy does not necessarily correlate with increased sensitivity and alignment quality (Yu et al., 2006;Sadreyev and Grishin, 2008). Therefore, to simultaneously improve all these qualities, we introduce simple adjustments to the predicted locationμ v (·) and scale parametersσ v (·) (v = A, B), as shown in Figure 3.
Observing a high correlation between the scale and the location parameters (Supplementary Section S6.3), we adjust the scale parameters as follows: The first equation eliminates the need to predictσ A . The second equation, although expressed non-linearly inσ B , employs one coefficient whose optimal value implies almost the same result as in the case of expressing the scale parameter linearly inμ B . We adjust the location parameters similarly: Note that the adjustmentsμ v andσ v (v = A, B) retain the dependence of the statistical parameters on the profile length and ENO and compositional similarity. They only scale and shift predictions made by the trained models. Conditional mean estimators (3) and (4) are used to combineμ v andσ v (v = A, B). The parameters a and b and the coefficients W = {gs, gi, gc, hs, hi, hc} in (5) and (6), which we refer to as the adjustment parameters, are optimized with respect to statistical accuracy and alignment quality and sensitivity (Supplementary Section S6.5). Here we note that the optimal balance is achieved with the compositional similarity between profiles having a large impact (a = b = 0.35) on conditional mean estimates.

Number of positive substitution scores
The number of positive profile-profile substitution scores in the alignment provides additional information to the alignment score. For example, the same alignment score can be the result of many weakly positive substitution scores or a few high scores. Therefore, it can be a useful indicator for estimating the significance of profile-profile alignments.
The distribution of the number of positive substitution scores can be accurately approximated by the negative binomial distribution (NBD) (Supplementary Section S6.4). However, aiming to reduce the overall complexity of the calculation of statistical significance, we propose a derived statistic ωn. It depends on the lengths of profiles being compared and its value (consequently, its significance) decreases as the alignment search space (the product of the profile lengths) increases: Here · denotes rounding to the nearest integer, l1 and l2 are the lengths of two profiles compared, ω is the number of positive substitution scores in the alignment, and c0 = 10 5 is the constant that is equal to the value of the denominator when l1 and l2 are slightly less than 50. Hence, the ωn statistic corresponds to the number of positive substitution scores when the search space approximately equals that of two profiles of length 50. The number of positive substitution scores becomes less informative as the search space increases, and, consequently, the value of the ωn statistic decreases. The distribution of ωn is approximately negative binomial (Supplementary Section S6.4). The complete procedure for statistical significance estimation is shown in Figure 3. p-values of alignment score and ωn, Pa and Po, are combined using the empirical Brown's method (Poole et al., 2016) (Supplementary Section S4.5).

E-value and its correction
The expected number of local alignments with a score greater than or equal to x, E = Kl1lN e −λx , increases as the size of the search space increases (Karlin and Altschul, 1990;Dembo et al., 1994;Waterman and Vingron, 1994b;Spang and Vingron, 2001). (E-value and p-value P are related by the equation P = 1−exp(−E).) Let EN be the E-value of an alignment with score x obtained from searching a database of size lN . Let also E N be the E-value of an alignment with the same score x (and profile composition) obtained by searching a database of size l N with the same query. Then, the relationship between EN and E N is (Yu et al., 2006) We use this relationship when compensating for a limited number of randomly generated profiles used in simulations. We refer to the results obtained by calculating the corrected E-value as the alternative model.

Results
This section presents results from the application of the proposed procedure (Figure 3) for estimating the statistical significance of profile-profile alignments. The NN models predicting statistical parameters were trained on the estimates obtained for profiles generated by Algorithm 1 with a noise level r = 0.03 and a fragment length s = 9, as determined in Section 2.4.

Statistical accuracy
We used COMER to search and align the profiles constructed for 5000 simulated Pfam (Finn et al., 2016) (v30.0) MSAs against the database of simulated profiles representing 4931 SCOPe (v2.03) domains. The profiles were randomized using Algorithm 1 to preserve length, ENO, and composition inherent in each of the Pfam MSAs and profiles constructed for the SCOPe domains (see Supplementary Section S5.2). Then, we calculated the fraction of queries with a p-value reported by COMER for their best match less than the specified p-value.  The results are shown in Figure 4. The new model for statistical significance estimation is more accurate than the model (Margelevičius and Venclovas, 2010) implemented in the previous version (Margelevičius, 2018) of the COMER method. The results also show that the statistics obtained from aligning profiles generated by randomly permuting MSA columns (s = 1, r = 0.01) lead to overestimation of statistical significance. This result can be accounted for by unrealistic representation of protein sequence families (profiles).
Finally, although the models with and without the ωn statistic taken into consideration achieve comparable statistical accuracy, the sensitivity and alignment quality obtained using the latter model (w/o ωn) are lower.

Sensitivity and alignment quality
We compared the performance of the new COMER version with that of its previous version (Margelevičius, 2018), where the new and previous versions differed only in how they estimated the statistical significance of alignments. For reference, we also provide results for two other profile-profile alignment methods, HHsearch (Söding, 2005) (v3.0.0) and FFAS (Jaroszewski et al., 2011). 4900 protein domains of the test dataset from the SCOPe database (v2.03) filtered to 20% sequence identity was used to evaluate performance. The test and training datasets shared no common folds.
Profiles were constructed using two categories of MSAs. The MSAs for each domain sequence were obtained by running PSI-BLAST (Altschul et al., 1997) (v2.2.28+) for six iterations and HHblits (Remmert et al., 2012) for three iterations, respectively, against a filtered UniProt database (Suzek et al., 2015). A pair of aligned domains (profiles) that belonged to the same SCOPe superfamily or shared statistically significant structural similarity (DALI (Holm et al., 2008) Z-score ≥ 2) was considered a true positive (TP). Aligned pairs that did not meet the above criteria but belonged to the same SCOPe fold were considered to have an unknown relationship and were ignored. Other aligned pairs were considered false positives (FPs). The sensitivity was also summarized using the ROCn score, which is the normalized area under the ROC curve up to n FPs. Alignment quality was evaluated by generating, using MODELLER (Šali and Blundell, 1993) (v9.4), protein structural models for each produced alignment. Statistically significant similarity between a model and the real structure, TM-score ≥ 0.4 (Zhang and Skolnick, 2004), was considered to correspond to a high-quality alignment (HQA). An alignment with a TM-score < 0.2, a characteristic value for a random pair, was assumed to be of low-quality (LQA).
Two evaluation modes were used to evaluate alignment quality. The local mode penalizes alignment overextension, whereas the global mode penalizes too short alignments (e.g., a few amino acids in length). These evaluation modes were also used to evaluate the quality of maximally extended alignments produced by COMER and HHsearch (option -mact set to 0). The evaluation of maximally extended alignments reveals how the quality of local alignments changes with their extension. It provides an indication of the quality of the match between the query and a database protein, which is also important in protein homology modeling.
(Details about the evaluation setting can be found in Supplementary Section S5.3.) Figure 5 and Table 1 reveal that the new statistical model lead to consistent improvement in both sensitivity and alignment quality and that most of the improvements are statistically significant.
Using the developed statistical model yielded an increase of up to 34.2% and 27.4% in the number of TPs and up to 43.9% and 61.8% in the number of HQAs. Table 2 shows that these percentage increases were achieved at a low false discovery rate (FDR). Hence, the largest relative improvements are expected for relationships detected at a high confidence level.
Table 2 also shows that the relative increases were greater when the profiles were constructed from the PSI-BLAST MSAs. A similar trend was also observed for the number of HQAs examined as a function of the number of alignments of inferior quality (IQAs) (Sup-  plementary Section S7.1). In this evaluation setting, the developed statistical model showed an increase of up to 102.8% and 193.3% in the number of HQAs.
MSAs obtained from a PSI-BLAST search contain more alignment errors than those built using HHblits. Larger relative improvements in performance achieved when using PSI-BLAST MSAs for profile construction, therefore, show a certain degree of robustness that the new model exhibit. We attribute this characteristic to the combination of statistical parameters dependent upon both profile length and ENO and compositional similarity. If a high alignment score of two unrelated profiles arises due to false positives in the corresponding MSAs, a high compositional similarity between the profiles counterbalances the statistical significance estimated solely based on profile length and ENO.

Application to pairwise profile HMM alignments
We applied the methodology for estimating statistical significance ( Figure 3) to pairwise profile HMM alignments produced by HHsearch. An improvement in both statistical accuracy and sensitivity and alignment quality (Supplementary Section S7.2) confirms the effectiveness of the developed methodology.

Example of reduced significance for a false positive
Reranking alignments using the proposed estimation of statistical significance has been shown to increase sensitivity and the rate of HQAs. This is demonstrated by an example.
Alpha-helical structure implies a high compositional similarity λu = 0.296 between the profiles constructed for the two domains (low values of λu correspond to high compositional similarity; Supplementary Sections S6.1 and S6.2). Significance estimation dependent upon compositional similarity allowed the new COMER version to correctly remove the alignment from the list of statistically significant alignments. In contrast, the alignment was considered significant by the previous COMER version. Note that both versions estimated the significance of the same alignment with the same score.

Discussion
Profiles represent sequence families, and this fact alone suggests that a profile contains information whose content largely depends on the family the profile describes. It means that the degree of (dis)similarity (distance) between unrelated sequence families strongly affects the distribution of scores of alignments between profiles that represent these families. The challenging aspects of characterizing the distribution of profile-profile alignment scores, however, are not limited to diversity across unrelated profiles. Similarity between SS, context and other predictions that accompany profiles complicate the characterization of alignment score distribution too.
The importance of a null model of profiles motivated us to develop an algorithm for generating random profiles. According to the algorithm, profiles are generated by concatenating fixed-length fragments sampled randomly from real unrelated profiles (seeds) with added noise. Controlling the noise level and fragment length allows to produce random profiles that are neither overly divergent nor share an excessive similarity. In this way, random profiles possess features characteristic to real profiles, but correlations that do not allow them to be considered unrelated do not dominate. Alignments between randomly generated profiles allowed us to determine the dependence of statistical parameters on profile length and ENO and also on compositional similarity between profiles.
The results suggest two important implications. First, profiles generated using long fragments (fragment length 9) represent real unrelated profiles much more accurately than do those obtained by randomly sampling profile columns. Statistics obtained from their alignments lead to higher statistical accuracy. On the other hand, randomly sampling of columns destroys higher-order dependencies inherent in real proteins and leads to the opposite result.
Second, the compositional similarity between profiles has a large impact on estimating the statistical significance of alignments. Improvements in sensitivity and alignment quality can be accounted for in part by that a high compositional similarity may indicate that the proteins share common structural elements or the profiles involve false positives. Significance estimation dependent upon compositional similarity, therefore, has a positive effect.
In fact, the issues of profile composition and random profile model are factors that hinder the application of techniques such as importance sampling (Park et al., 2009) to accurately estimating statistical parameters. While composition can be measured in different ways, choosing a random profile model is more complicated because every profile represents a model. For example, the composition of two profiles, one of them obtained by randomly rearranging the positions of the other, will be the same. However, the distribution of alignment scores of profiles generated using the profile with rearranged positions as a seed (model) will differ from that obtained for profiles generated using the other profile as a seed. In this study, close match to the distribution of alignment scores of real unrelated reference profiles constituted one of the criteria for random profile model selection. However, a supplementary approach might be beneficial.
Based on the above considerations and the results obtained by generating random profiles using one seed profile, we suspect that further improvements may result from introducing a new measure of profile "sequence affinity" (as opposed to the quantitative measure of the effective number of observations). The benefit may come from including into the model the statistical parameters dependent upon the new measure through, e.g., the application of the conditional mean estimator. The present study already demonstrated the conditional mean estimator to be effective when the calculation of compositional similarity or divergence between profiles was in use. And it provides possibilities for further improvements.