# R: statistical test to identify samples with too high variability

I would like to develop a test to identify which variables in my data set have a variation higher than the "average variability".

I'm struggling with that since days, and I also tried in vain to look for help in other forums.

I have data from biological experiments, that look like this:

v1 2 1.8 1.5 1.9 2.1 1.78 1.95 2.0 2.1
v2 2 100 -5.2
v3 1 -1.3 -2 2.3
v4 1 1.5 1.6 1.9 2.1 2.0 2.4 -1.1 2.3 1.5 1.6 1.9 1.8 1.6


These represent gene expressions. Now, I would expect that all values of each variable(genes) are more or less similar, since the values are repeat measurements of the same gene.

Having a variable with such a huge difference, as v2 , doesn't have sense, because the repeated measurements should give consistent values. Therefore, it has to come from a methodological error and the variable (gene) has to be discarded.

I was looking for a method (possible a statistical test) in R which could identify the "average variability" among my samples and report me which variables (genes) have a variability significantly greater. This means that for these genes my data are not good enough to estimate the expression, and I have to discard them.

I would really appreciate any suggestion/links/advice/methods on test I could use for my purpose.

• ok, thank you and sorry, but also there I did not manage to find much help. – efrem Jul 3 '14 at 10:42
• If you think it is the case, feel free to put this off. – efrem Jul 3 '14 at 10:50

I just looked at this.

My approach was:

• compute the mean, standard deviation, and count for each set of samples
• compute the critical t-threshold given alpha, the sample size, and the nature of the fit (quadratic). I was using excel so I used "T.inv".
• transform the data by subtracting the mean, then dividing by the standard deviation, then comparing the absolute value to the t-threshold.
• If it is above the threshold then it is classified as an outlier

Note: alpha is a parameter. If you want to make your fit "wider" then use a smaller value. If you want more data to be classified as possible outlier then use a higher value. It is exceptionally good if you can take the time to understand what "alpha" means in the statistical sense of this threshold.

I notice you have rows with 3 samples - that is dangerous:

Having two samples and computing the standard deviation is like having one sample and computing the mean. The math gives you a number, but it is as sample-sparse as mathematics can go and still give a value - it is on the edge of the cliff of oblivion and is not very informative. Get more samples.

There are rules of thumb that say 5, 10, 30, 100 or 300 are sufficient. If you are going below 5 then you had best have a great defense for why the math isn't bad.

The "average variability" that you want to measure, should translate in Standard Deviation for statistics. It's pretty easy to compute STD in R, so look up the definition of Standard Deviation on google to see if it matches with what you want to find.

• I referred to the "average variability" as a value which estimate how much each samples can vary modeled from the variability/dispersion of the single variables. I do not refer to the variability of each variable. – efrem Jul 3 '14 at 10:40
• Yes, you are right. But we are still talking about different things. I want to create a model considering std.v1, std.v2, std.v3..., which I can use to predict that with a certain confidence a variable in this dataset has a std in a specific range. Samples with high std will be out of the model, and I want to identify and remove them. – efrem Jul 3 '14 at 10:49
• Compute std of each sample. Remove samples with high stds that mean there was some error. Computation of std is easy in R. You'll have to define somehow the accepted std, which is specific to your application. Your problem is the programmatistical implementation of 1.find std of each sample, 2. remove rotten samples. ? – a_kats Jul 3 '14 at 11:03
• @efrem: To translate "a_kats"'s suggestion in R code: if you have all your "v"s in a list, you could sapply(my_list_of_v, var) and then call boxplot.stats with a suitable for your needs "coef" argument. E.g. boxplot.stats(sapply(my_list_of_v, var), coef = 1.5)$out returns no outliers while boxplot.stats(sapply(my_list_of_v, var), coef = 0.9)$out returns. – alexis_laz Jul 3 '14 at 11:20
• thank you. I guess that the boxplot calculate outliers based on quantiles and assuming normally distributed data. For me is not always the case, so I was looking for something which could consider the trend of the data as well. – efrem Jul 3 '14 at 11:53