You may want to consider using an intraclass correlation (ICC) for this purpose. At its heart, this is an analysis of variance, like an ANOVA. Basically, the idea is to see what proportion of the total variance can be attributed to the CupCake ID. To the extent that there are other sources of variance in your data (such as differences between raters), the ICC will drop. It can also be interpreted as the correlation among ratings of the same cupcake --- if your raters always perfectly agreed, this would be 1.
The psych package in R provides a nice ICC function, and excellent help documentation that briefly describes not only the function, but a more general explanation of how ICCs work and what they mean.
An example in R
First, here's some made-up data matching the structure you describe. I've created a set of "true" ratings for each cupcake. The ratings from each rater mostly match the true ratings, but for each rater, I replaced a random 1/3 of their ratings with a new random rating, so that there is some disagreement.
set.seed(24601)
true_ratings <- sample(1:4, 300, replace = TRUE)
Rater1 <- true_ratings
Rater1[sample(1:300, 100)] <- sample(1:4, 100, replace = TRUE) # replace 1/3 with random ratings
Rater2 <- true_ratings
Rater2[sample(1:300, 100)] <- sample(1:4, 100, replace = TRUE) # replace 1/3 with random ratings
Rater3 <- true_ratings
Rater3[sample(1:300, 100)] <- sample(1:4, 100, replace = TRUE) # replace 1/3 with random ratings
Rater4 <- true_ratings
Rater4[sample(1:300, 100)] <- sample(1:4, 100, replace = TRUE) # replace 1/3 with random ratings
Rater5 <- true_ratings
Rater5[sample(1:300, 100)] <- sample(1:4, 100, replace = TRUE) # replace 1/3 with random ratings
data <- data.frame(MachineID = rep(1:100, each=3),
CupCakeID = 1:300,
Rater1 = Rater1,
Rater2 = Rater2,
Rater3 = Rater3,
Rater4 = Rater4,
Rater5 = Rater5)
Here's what the first 6 lines of that look like:
> head(data)
MachineID CupCakeID Rater1 Rater2 Rater3 Rater4 Rater5
1 1 1 2 2 2 2 3
2 1 2 1 1 1 1 3
3 1 3 4 3 3 3 3
4 2 4 4 4 4 4 4
5 2 5 1 3 2 3 3
6 2 6 2 3 4 2 2
I am only concerned about the objectivity of the tester and not the
quality of the cupcake machines.
If you really don't care about the machine and you're confident the machines are all interchangeable, you can discard it from the model completely and just look at cupcakes and raters.
library(psych)
ICC(data[, 3:7]) # just using the rater columns, each row is one cupcake
This provides the output for 6 different versions of the ICC.
Call: ICC(x = data[, 3:7])
Intraclass correlation coefficients
type ICC F df1 df2 p lower bound upper bound
Single_raters_absolute ICC1 0.43 4.8 299 1200 0 0.38 0.49
Single_random_raters ICC2 0.43 4.8 299 1196 0 0.38 0.49
Single_fixed_raters ICC3 0.43 4.8 299 1196 0 0.38 0.49
Average_raters_absolute ICC1k 0.79 4.8 299 1200 0 0.75 0.83
Average_random_raters ICC2k 0.79 4.8 299 1196 0 0.75 0.83
Average_fixed_raters ICC3k 0.79 4.8 299 1196 0 0.75 0.83
Number of subjects = 300 Number of Judges = 5
The correct ICC to interpret for your situation depends on how you're thinking about the data generating process. If you want to consider your raters as a random effect, then the ICC2 (either Single_random_raters or Average_random_raters) makes the most sense in your case.
The ICC function is handy, but it really is just a ratio of variances, so it's not hard to calculate by hand if you want to allow for a different model specification. For example, you may decide you want to allow for the possibility of variance in ratings due to machine. To flexibly calculate ICC for whatever model structure you need, use the variance estimates from a mixed effects model. Here is an accessible discussion of ICCs calculated from mixed effects models, a useful answer here on SE, and several examples of ICCs calculated from mixed effects models.
