BSig: evaluating the statistical significance of biclustering solutions

Abstract

Statistical evaluation of biclustering solutions is essential to guarantee the absence of spurious relations and to validate the high number of scientific statements inferred from unsupervised data analysis without a proper statistical ground. Most biclustering methods rely on merit functions to discover biclusters with specific homogeneity criteria. However, strong homogeneity does not guarantee the statistical significance of biclustering solutions. Furthermore, although some biclustering methods test the statistical significance of specific types of biclusters, there are no methods to assess the significance of flexible biclustering models. This work proposes a method to evaluate the statistical significance of biclustering solutions. It integrates state-of-the-art statistical views on the significance of local patterns and extends them with new principles to assess the significance of biclusters with additive, multiplicative, symmetric, order-preserving and plaid coherencies. The proposed statistical tests provide the unprecedented possibility to minimize the number of false positive biclusters without incurring on false negatives, and to compare state-of-the-art biclustering algorithms according to the statistical significance of their outputs. Results on synthetic and real data support the soundness and relevance of the proposed contributions, and stress the need to combine significance and homogeneity criteria to guide the search for biclusters.

Appendices

Appendix 1: Relevance of biclustering with flexible coherency

High-dimensional biomedical and social data is characterized by the presence of biclusters with flexible coherency assumptions (Table 6). Table 6 motivates the relevance of such biclusters, highlighting data contexts where their discovery is relevant for different purposes.

Table 6 Relevance of non-constant biclusters for biomedical and social data analysis

1.1 Appendix 2: Continuous adjustment factors

The proposed statistical tests can be extended to support biclusters with continuous adjustment factors. Consider the $\{1.3,2.2,1.7\}$ combination of values, a continuous range of coherent values under additive or multiplicative assumptions can be generated based on the exploration of $\gamma$ factors (e.g. shifting $\gamma \in [-1.3,1.8]$ or scaling $\gamma \in [0,1.8]$ factors for values $a_{ij}\in [0,4]$). In order to robustly compute the $p_{\varphi _B}$ probability of additive and multiplicative models, and subsequently of symmetric and plaid models (with underlying additive/multiplicative assumptions), we propose a technique based based on the integral of the product of (either slided or scaled) density probability functions.

1.2 Continuous ranges of shifting factors

Consider the additive coherency assumption. Let the maximum and minimum observed values for a particular row $x_i\in I$ of a bicluster to be, respectively, $max_{J|x_i}$ and $min_{J|x_i}$. Also consider the range of real-values of the matrix A to be $[min_{A},max_{A}]$. Then for a particular pattern $\mathbf {J}|\mathbf {x}_i$ the shifting factors are defined by the interval $\gamma \in [\gamma _1=-(min_{A}-min_{J|x_i}),\gamma _2=max_{A}-max_{J|x_i}]$. The probability of a particular value $a_{ij}$ to occur under this shifting interval is:

$$\begin{aligned} \int _{a_{ij}+\gamma _1}^{a_{ij}+\gamma _2}f(x)=\int _{\gamma _1}^{\gamma _2}f(x+a_{ij}) \end{aligned}$$

(A1)

where f(x) is the distribution function that approximates $a_{ij}$ values. This calculus assumes that the range of observed values $\hat{A}$ are linearly adjusted to guarantee an unitary coherency strength $\delta \approx 1$. The probability of two values $a_{ij}$ ($a_1$) and $a_{i(j+1)}$ ($a_2$) to occur under this shifting interval is not simply the product of their individual probabilities since a simple product would allow for non-coherent values (e.g. $\{a_{1}+\gamma _1,a_{2}+\gamma _2/2\}$). In order to correctly account for the combination of values with continuous shifting ranges, the distribution functions need to be aligned by the target column value and multiplied. The resulting function delivers the product of the individual probabilities. Finally, the area behind this curve between $\gamma _1$ and $\gamma _2$ values is computed in order to retrieve a estimate of the probability $p_{\varphi _B}$ for the Binomial tail calculus. This strategy is illustrated in Fig. 14, under the assumption that the values in A are either described by an Uniform or Gaussian distribution. Given $\varphi _{B}^i=\{a_{i1},..,a_{im}\}$ combination of values, $p_{\varphi _{B}^i}$ can be approximated by:

$$\begin{aligned} \int _{\gamma _1}^{\gamma _2}{\varPi }_{j=1}^m f(x+a_{ij}) \end{aligned}$$

(A2)

In order to compute this probability efficiently we propose the calculus of its approximate area by interpolating 100 points between $\gamma _1$ and $\gamma _2$.

Fig. 14

Illustrative integral of the product of slided density functions to assess biclusters with continuous ranges of shifting factors

1.3 Continuous ranges of scaling factors

The probability of occurrence of a combination of real values $\varphi _{B}^i$ on the $i\hbox {th}$ row of a bicluster under a multiplicative coherency across rows can be approximated using similar principles to the ones proposed in previous section. Considering $max_i$ and $min_i$ to be the maximum and minimum values of a given row $x_i$ and $\bar{A}$ to be the range real values in A. When only positive values are allowed, the scaling range is $[\gamma _1=0,\gamma _2=\bar{A}/max_i]$. When negatives values are allowed the scaling range is given by $[\gamma _1=-d,\gamma _2=d]$ where $d=max(max_i,-min_i)/(\bar{A}/2)$.

The probability of multiple values to occur is given by the integral of the product of the size-adjusted density functions for the $[\gamma _1,\gamma _2]$ interval. Why the size adjustment is necessary? Consider the pair of observed values $\{a_1=1,a_2=2.5\}$ and the scaling range to be $\gamma \in [0,1]$. This means that the density function to estimate the $a_1$ value is considered for the interval [0,1], while the density function to estimate $a_2$ is considered over [0,2.5]. Therefore, the density functions need to be normalized with regards to their size: $f(x/a_1)$ and $f(x/a_2)$. Given $\varphi _{B}^i=\{a_{i1},..,a_{im}\}$ combination of values for i row, $p_{\varphi _{B}^i}$ can be approximated by:

$$\begin{aligned} \int _{c_1}^{c_2}{\varPi }_{i=1}^n f(x/a_i) \end{aligned}$$

(A3)

Similarly, an efficient computation of the (A3) integral calculus is made available recurring to interpolation whenever the multiplication of the inputted density functions is complex. This strategy is illustrated in Fig. 15, under the assumption that the values of the A matrix are either described by a single Uniform or Gaussian distribution.

Fig. 15

Illustrative integral of the product of scaled density functions to assess biclusters with continuous ranges of scaling factors

Table 7 Impact of data size and dimensionality on the expected minimum number of observations $(\mathrm{\mathbf{n}}_{min})$ in biclusters with continuous adjustment factors to guarantee their statistical significance (assuming a $\delta =0.2$ coherency strength, uniform background values, and additive and multiplicative coherencies with varying ranges of allowed shifts/scales). Algorithm 1 was applied to compute statistical significance

Appendix 3: Complementary results

Table 7 shows how the required minimum of rows that guarantee the statistical significance of a real-valued bicluster with continuous shifts/scales varies with the number rows and columns of the input dataset. Two major observations can be retrieved. Both the size and dimensionality of data affect the significance levels, being the effect of varying the size of data clearly more accentuated since the assessment was applied over biclusters with coherency across rows. Second, the observed pattern also largely determines the computed significance levels as it determines the range of allowed shifts and scales (Table 8). Understandably, the larger the allowed range, the higher is the probability of a bicluster pattern to occur and thus the higher is the number of minimum rows in the bicluster to guarantee its significance.

Table 8 Expected probability of different patterns to occur in biclusters with continuous shifts and scales from data with approximately uniform distribution of values ($a_{ij}\in$[0,1])

Fig. 16

Impact of coherency strength, pattern length, data size/dimensionality on the expected minimum number of bicluster’s rows that guarantee its statistical significance

Figure 16 provides the graphical representation of the results gatheres throughout Tables 1, 2 and 3, thus showing the expected minimum number of rows in a bicluster that guarantees its significance for varying: coherency assumption, pattern expectations $\varphi _B$, coherency strength $|\mathcal {L}|$, pattern length m, data size N, and data dimensionality M.