1
$\begingroup$

I have two algorithms running on 100 benchmark instances.
I want to test if there are significant differences in their mean average error.

A significant part of those two algorithms is not deterministic.

My data with only one run looks like this:

+--------------------+-------+-----------+
| benchmark instance | error | algorithm |
+--------------------+-------+-----------+
| BENCH1             | 1.1   | ALGO1     |
+--------------------+-------+-----------+
| BENCH2             | 5.4   | ALGO1     |
+--------------------+-------+-----------+
| BENCH3             | 7.8   | ALGO1     |
+--------------------+-------+-----------+
| ...                | ...   | ...       |
+--------------------+-------+-----------+
| BENCH1             | 2.5   | ALGO2     |
+--------------------+-------+-----------+
| BENCH2             | 4.4   | ALGO2     |
+--------------------+-------+-----------+
| BENCH3             | 7.7   | ALGO2     |
+--------------------+-------+-----------+
| ...                | ...   | ...       |
+--------------------+-------+-----------+
| BENCH1             | 1.9   | ALGO1     |
+--------------------+-------+-----------+
| BENCH2             | 5.2   | ALGO1     |
+--------------------+-------+-----------+
| BENCH3             | 6.1   | ALGO1     |
+--------------------+-------+-----------+

In this case I use 2-sample t-test.

t.test(error ~ algorithm, data = combined.data)

If I run the 100 tests for a second or a third time for each algorithm, do I need to use paired t-test or not?

t.test(error ~ algorithm, data = combined.data, paired = TRUE)

I use as a rule of thumb that when measurements on same group of objects are made twice, then in order to compare the difference paired T-test is used. Is that correct in that case? Should a repeated test on the same instance with the same algorithm considered as the same object?

The fact is that I need to run each algorithm many times on the same instances in order to be more precise. Some facts like CPU lagging or another app use some portion of the CPU or the randomized nature of some commands inside the algorithm may give an error of lets say 11% and the second time the same algorithm on the same instance will give 9%. If I run it 50 times I can say for sure that the quality of this algorithm ranges from 9-11%.

The final goal is to find out which algorithm performed better by comparing the error column.

Thanks

--------------------------------------------------------

Addition,

Data with 5 repetitions * 7 instances * 2 Algorithms.

+----------+-------+---------+
| bench_id | Algo  |  Score  |
+----------+-------+---------+
| BENCH25  | ALGO1 |  4.7601 |
| BENCH14  | ALGO1 |  2.9453 |
| BENCH13  | ALGO1 |  3.0815 |
| BENCH10  | ALGO1 |  6.3272 |
| BENCH6   | ALGO1 |  8.7662 |
| BENCH101 | ALGO1 |  6.7073 |
| BENCH200 | ALGO1 |  7.7389 |
| BENCH25  | ALGO1 |  4.7404 |
| BENCH14  | ALGO1 |  2.8845 |
| BENCH13  | ALGO1 |  4.0013 |
| BENCH10  | ALGO1 |  6.3065 |
| BENCH6   | ALGO1 |  9.6564 |
| BENCH101 | ALGO1 |  7.1119 |
| BENCH200 | ALGO1 |   8.539 |
| BENCH25  | ALGO1 |  4.4898 |
| BENCH14  | ALGO1 |  2.8129 |
| BENCH13  | ALGO1 |  3.7915 |
| BENCH10  | ALGO1 |  7.1798 |
| BENCH6   | ALGO1 |  8.1818 |
| BENCH101 | ALGO1 |  5.8875 |
| BENCH200 | ALGO1 |     8.6 |
| BENCH25  | ALGO1 |  3.0398 |
| BENCH14  | ALGO1 |  2.8487 |
| BENCH13  | ALGO1 |  4.4973 |
| BENCH10  | ALGO1 |  7.6505 |
| BENCH6   | ALGO1 |  7.9472 |
| BENCH101 | ALGO1 |  6.4629 |
| BENCH200 | ALGO1 |  9.0849 |
| BENCH25  | ALGO1 |  4.3706 |
| BENCH14  | ALGO1 |  2.7197 |
| BENCH13  | ALGO1 |  4.1304 |
| BENCH10  | ALGO1 |  6.4852 |
| BENCH6   | ALGO1 |  8.6536 |
| BENCH101 | ALGO1 |  6.7073 |
| BENCH200 | ALGO1 |  9.2053 |
| BENCH25  | ALGO2 |  5.1921 |
| BENCH14  | ALGO2 |  4.4924 |
| BENCH13  | ALGO2 |  6.1262 |
| BENCH10  | ALGO2 |  5.9123 |
| BENCH6   | ALGO2 |  9.0649 |
| BENCH101 | ALGO2 |  6.9505 |
| BENCH200 | ALGO2 | 10.0394 |
| BENCH25  | ALGO2 |   6.054 |
| BENCH14  | ALGO2 |  4.4786 |
| BENCH13  | ALGO2 |   5.641 |
| BENCH10  | ALGO2 |  9.8324 |
| BENCH6   | ALGO2 |  9.9433 |
| BENCH101 | ALGO2 |  7.6724 |
| BENCH200 | ALGO2 |  8.8431 |
| BENCH25  | ALGO2 |  6.3283 |
| BENCH14  | ALGO2 |  3.9496 |
| BENCH13  | ALGO2 |  6.8156 |
| BENCH10  | ALGO2 |  6.8474 |
| BENCH6   | ALGO2 |  7.9981 |
| BENCH101 | ALGO2 |  8.3048 |
| BENCH200 | ALGO2 | 10.5675 |
| BENCH25  | ALGO2 |  5.2637 |
| BENCH14  | ALGO2 |  4.4197 |
| BENCH13  | ALGO2 |  5.5159 |
| BENCH10  | ALGO2 |   8.963 |
| BENCH6   | ALGO2 |    8.39 |
| BENCH101 | ALGO2 |  7.5928 |
| BENCH200 | ALGO2 |  9.9212 |
| BENCH25  | ALGO2 |  5.4682 |
| BENCH14  | ALGO2 |  5.0907 |
| BENCH13  | ALGO2 |  5.5636 |
| BENCH10  | ALGO2 |  8.8328 |
| BENCH6   | ALGO2 |  8.2008 |
| BENCH101 | ALGO2 |  7.2727 |
| BENCH200 | ALGO2 |  9.5646 |
+----------+-------+---------+

Using 1st approach like this:

b_var <- 1/var(combined.data$score)
res <- lmer(
  score ~ algo + (1 | bench_id), combined.data,
  weights = 1 / b_var)
summary(res)

I get Error variable lengths differ (found for '(weights)')


Using 2nd approach I get:

res1 <- lmer(
  score ~ algo + (1 | bench_id), combined.data)
summary(res1)
Linear mixed model fit by REML ['lmerMod']
Formula: score ~ algo + (1 | bench_id)
   Data: combined.data

REML criterion at convergence: 183.2

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.91976 -0.48884  0.00928  0.59819  2.47335 

Random effects:
 Groups   Name        Variance Std.Dev.
 bench_id (Intercept) 4.3060   2.0751  
 Residual             0.5283   0.7269  
Number of obs: 70, groups:  bench_id, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept)   5.9518     0.7939   7.497
algoALGO2     1.2229     0.1738   7.038

Correlation of Fixed Effects:
            (Intr)
algoALGO2   -0.109


Using 3rd approach I get:

res2 <- lmer(
  score ~ algo + (0 + factor(algo) | bench_id), combined.data)
summary(res2)
Linear mixed model fit by REML ['lmerMod']
Formula: score ~ algo + (0 + factor(algo) | bench_id)
   Data: combined.data

REML criterion at convergence: 177.5

Scaled residuals: 
     Min       1Q   Median       3Q      Max 
-2.91175 -0.51624 -0.00599  0.60168  2.73685 

Random effects:
 Groups   Name                Variance Std.Dev. Corr
 bench_id factor(algo)ALGO1  5.2709   2.296        
          factor(algo)ALGO2 3.4470   1.857    1.00
 Residual                     0.4816   0.694        
Number of obs: 70, groups:  bench_id, 7

Fixed effects:
            Estimate Std. Error t value
(Intercept)   5.9518     0.8756   6.797
algoALGO2     1.2229     0.2347   5.210

Correlation of Fixed Effects:
            (Intr)
algoALGO2   -0.768
convergence code: 0
singular fit


Using t-test (paired):

data:  score by algo
t = -7.0242, df = 34, p-value = 4.165e-08
alternative hypothesis: true difference in means is not equal to 0
95 percent confidence interval:
 -1.5766716 -0.8690655
sample estimates:
mean of the differences 
              -1.222869


Using Wilcoxon test:

wilcox.test(score ~ algo, data = combined.data)

  Wilcoxon rank sum test with continuity correction

data:  score by algo
W = 425, p-value = 0.02805
alternative hypothesis: true location shift is not equal to 0
$\endgroup$
  • $\begingroup$ You have two samples, algo1 and algo2. But they are linked by the benches assuming bench1 is the same task for both algorithms. So you should use a paired t-test. If the benches are arbitrary, then independent samples. $\endgroup$ – Heteroskedastic Jim Nov 15 '18 at 12:56
  • $\begingroup$ @HeteroskedasticJim would it have any negative impact if I run the test 2-3 times and add those rows? I mean would it contribute something or make it worse? $\endgroup$ – Antonis Nov 15 '18 at 12:58
  • $\begingroup$ Why would you do this? Help explain. $\endgroup$ – Heteroskedastic Jim Nov 15 '18 at 12:59
  • 1
    $\begingroup$ Because when I run only one test on one instance with Algorithm A for example. some facts like CPU lagging or another app use some portion of the CPU or the randomized nature of some commands inside the algorithm may give an error of lets say 11% and the second time the same algorithm on the same instance will give 9%. If I run it 50 times I can say for sure that the quality of this algorithm ranges from 9-11% $\endgroup$ – Antonis Nov 15 '18 at 13:04
  • 1
    $\begingroup$ You should modify your question to make this part clear. The most natural approach would be a mixed effects regression. I'll write an answer. $\endgroup$ – Heteroskedastic Jim Nov 15 '18 at 13:07
2
$\begingroup$

The most natural approach to handle your situation would be a mixed effects or multilevel regression. You can also perform an inverse variance weighted regression.

You have 100 benchmarks for both algorithms. Additionally, you plan to attempt each benchmark about 50 times to account for randomness in the benchmark results.

Inverse variance weighting

I'll begin with the inverse variance approach since it is easier. Assuming that the data for the 50 instances of a benchmark for an algorithm are reasonably normal, then the data can be summarized by its mean and variance. So for each benchmark and algorithm, calculate the mean and the variance of the 50 instances. And you will arrive at 100 average benchmarks for each algorithm with 100 variances.

The averages become the outcome for a regression model where you correlate the scores by benchmark to account for the fact that each benchmark produced two averages. You also need to supply the variances as weights to the regression. Assume the benchmark averages are score, benchmark unique identifier is bench_id, algorithm indicator or dummy variable is algo, and benchmark variances are b_var, all contained in the dataset, dat. dat should be 200 rows long, then:

library(lme4)
res <- lmer(
  score ~ algo + (1 | bench_id), dat,
  weights = 1 / b_var)
summary(res)

Should work. You'll get the difference between algorithms from its coefficients, and you'll get a bench random intercept standard deviation which I'll call $v$, and a residual standard deviation which I'll call $s$. If you perform:

$$\frac{v^2}{v^2+s^2}$$

You'll get the average correlation between measurements from the same benchmark, this is known as the (residual) ICC. The syntax should be:

as.numeric(VarCorr(res)) / (as.numeric(VarCorr(res)) + sigma(res) ^ 2)

Complete multilevel setup

This time, you have all the data in long form. Each row pertains to one of 50 measurements for an algorithm for one of 100 benchmarks. So your number of rows should be $10000\ (50\times 100\times 2)$. The model you'll run still has a similar syntax to the 10000 long dataset:

res1 <- lmer(
  score ~ algo + (1 | bench_id), dat.longer)
summary(res1)

I would prefer to run:

res2 <- lmer(
  score ~ algo + (0 + factor(algo) | bench_id), dat.longer)
summary(res2)

This would allow you to observe the correlation between the scores from both algorithms and test whether the algorithms vary differently.


Of the three models, I prefer the first and the last. The first is a simple approach and if the normality assumption for the 50 instances per benchmark per algorithm is satisfied, then it is reasonable and defensible. Alternatively, the final model uses all the data in the most natural way. However, it assumes homogeneity of variance which the first approach does not. Final model assumes the variance of the 5000 measurements within each algorithm is the same across the 100 benchmarks.

Hope this helps. It might be useful to pick up a basic guide to multilevel modeling using lme4 if you are not familiar with such models. For your problem, especially if you stick with the first approach, an introductory guide is sufficient.


EDIT

Assuming the dataset is called data.long, then:

# to create dataset of mean scores
dat <- aggregate(score ~ algo + bench_id, data.long, mean)
dat$b_var <- aggregate(score ~ algo + bench_id, data.long, var)$score

Then use dat for model 1 syntax above. This dataset will have only 14 rows since you currently have only 7 bench_ids.

Try model 3 when you get more data with your combined.data.

$\endgroup$
  • $\begingroup$ Thank you so much for your thorough answer. I need to wait to accept it because I'm a new user. $\endgroup$ – Antonis Nov 15 '18 at 14:40
  • $\begingroup$ @Antonis glad it helps. Welcome to the site. I'm used to new users not getting the formatting correct so I assumed you'd been here longer. $\endgroup$ – Heteroskedastic Jim Nov 15 '18 at 14:42
  • $\begingroup$ I just wanted to give some feedback and ask some guidance. I hope I'm not bothering you. I run 7 instances for five times for both algos and I got back 70 lines. I tried all the models. Sadly only the second (the one that you do not recommend worked). I must do something wrong. At the first model I calculate weights like this: 1 / var(mydata$scores) and when I run the code it returns variable lengths differ (found for '(weights)'). When I run the 3rd model it returns singular fit. Does this imply that the algorithms do not differ? $\endgroup$ – Antonis Nov 15 '18 at 20:55
  • $\begingroup$ @Antonis see the edit at the bottom of the question. $\endgroup$ – Heteroskedastic Jim Nov 16 '18 at 3:21

Your Answer

By clicking "Post Your Answer", you acknowledge that you have read our updated terms of service, privacy policy and cookie policy, and that your continued use of the website is subject to these policies.

Not the answer you're looking for? Browse other questions tagged or ask your own question.