Suppose you ordered 50 different chemical products from a vendor.
The products are delivered to you in vials, each labelled with a unique barcode, but completely identical in all other respects.
The vendor sends you a (text) file where each barcode is matched to each product identifier.
The identifiers of the products you ordered are (I'm using R):
p <- paste("P",1:50,sep="_")
Each product has a property $M$, which is known with certainty.
You have your own table matching each product identifier to its $M$.
So, after using the vendor's file to match each barcode to each identifier, you match the identifier to $M$, and get the following:
pM <- c(136L, 137L, 107L, 84L, 89L, 86L, 96L, 143L, 100L, 93L, 84L, 115L, 117L, 89L, 144L, 87L, 148L, 84L, 121L, 96L, 149L, 130L, 115L, 122L, 135L, 128L, 145L, 88L, 133L, 109L, 102L, 108L, 119L, 87L, 126L, 148L, 89L, 87L, 115L, 142L, 118L, 92L, 96L, 118L, 81L, 148L, 121L, 95L, 140L, 118L)
You can measure $M$ very accurately, although it is expensive.
So, to test the correctness of the delivery without having to test all the products, you take a random sample from, say the vials bearing the barcodes that were matched to products P_4, P_15, P_43, and measure $M$ for them.
You expect to find:
pM[c(4,15,43)]
#[1] 84 144 96
but instead of (84, 144, 96), you find (148, 144, 96), i.e. the first sample is wrong.
Assuming that all the ordered products were delivered (i.e. that if all the samples were measured, all the expected $M$'s would be found, but not necessarily in the expected vials), the only possibility is that the data were scrambled, so some barcodes were wrongly matched to the products.
The question is: given your observation of 1 wrong sample out of 3 tested, what inferences can you make on the possible scrambling of the data, i.e. on the most likely number of samples that were swapped, if any?
I could not immediately place this problem theoretically, so I decided to try a simulation.
Essentially, for each n from 1 to 50 I randomly sampled 10000 times n elements of pM and swapped them. In each of these cases, I checked whether the actually observed situation (sample 4 showing the wrong $M$, samples 15 and 43 showing the correct $M$) was reproduced, and I counted how many times this happened over the 10000 repeats. Then I plotted these counts vs n.
probdist_dens = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; matches = (pMsw[c(4,15,43)] == pM[c(4,15,43)]); prod(matches==c(FALSE,TRUE,TRUE))})); probdist_dens = rbind(probdist_dens,data.frame(n=n,dens=dens)) }
plot(probdist_dens,type="l")
I don't know what to make of this, or even if it makes sense at all.
Can I conclude that the most likely number of swapped samples is close to 20?
The expected value seems to be close to 21 (which in this case means 20 because of the way I swapped the samples, just reversing their positions):
with(probdist_dens,sum(n*dens)/sum(dens))
#[1] 20.79004
If I test the case of finding the expected $M$ for all 3 samples, I get:
probdist_dens2 = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; matches = (pMsw[c(4,15,43)] == pM[c(4,15,43)]); prod(matches)})); probdist_dens2 = rbind(probdist_dens2,data.frame(n=n,dens=dens)) }
plot(probdist_dens2,type="l")
which seems to imply that even when all 3 samples seem correct, there could still be scrambling (the expected value is around 10).
And the case of not finding the expected $M$ for any of the 3 samples:
probdist_dens3 = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; matches = (pMsw[c(4,15,43)] == pM[c(4,15,43)]); prod(!matches)})); probdist_dens3 = rbind(probdist_dens3,data.frame(n=n,dens=dens)) }
plot(probdist_dens3,type="l")
(the expected value is around 41).
Maybe one reassuring fact is that when I expand the number of tested samples, but still maintaining the same proportion of mismatches, e.g. finding 2 incorrect samples on testing 6:
probdist_dens4 = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; matches = (pMsw[c(4,15,20,35,43,49)] == pM[c(4,15,20,35,43,49)]); prod(matches==c(FALSE,TRUE,FALSE,TRUE,TRUE,TRUE))})); probdist_dens4 = rbind(probdist_dens4,data.frame(n=n,dens=dens)) }
plot(probdist_dens4,type="l")
the expected value is about in the same place, but the distribution is more narrow.
Similarly, for the case where all samples have the incorrect $M$, if I take the number of tested samples to 6:
probdist_dens5 = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; matches = (pMsw[c(4,15,20,35,43,49)] == pM[c(4,15,20,35,43,49)]); prod(!matches)})); probdist_dens5 = rbind(probdist_dens5,data.frame(n=n,dens=dens)) }
plot(probdist_dens5,type="l")
with(probdist_dens5,sum(n*dens)/sum(dens))
#[1] 44.59095
Do you think my simulations address the problem I'm trying to solve?
(which BTW is a real situation I am tackling at work at the moment, not a thought experiment; only with many more chemical products, about 2000).
Thanks!
EDIT (following up from whuber's feedback)
I think whuber's main point is: the fact that in the sample of 3 vials that were measured, the one labelled P_4 (expected $M$ = 84) has instead $M$ = 148, tells us something about the remaining (non-sampled) vials, too.
There are 3 vials with $M$ = 84, whose correct labels should be:
p[pM==84]
#[1] "P_4" "P_11" "P_18"
We know that label P_4 was incorrectly placed on a vial containing a chemical product with different $M$ (184).
This implies that, in the set of vials that we did not sample, there is a vial containing product P_4, i.e. with $M$ = 84, but bearing the wrong label.
The wrong label could in theory be P_11 or P_18, i.e. of a product with the same $M$ as P_4, but that does not matter, because it would still imply that either P_11 or P_18 have in turn a label of a product with a non-matching $M$.
So the finding of 1 $M$ mismatch in the 3 tested samples implies by necessity the existence of another $M$ mismatch in the non-tested samples (a sample with $M$ 84 but labelled with a product identifier corresponding to a different $M$).
I am not entirely sure how to use this in practice, mainly because for the simulation (at least in the way I built it) one has to fix the set of tested samples, and as there are duplicated $M$'s, it is not indifferent which other mismatched sample is chosen. Maybe the other mismatched sample can be randomly sampled each time among the suitable complementary set.
probdist_dens_new = data.frame(); for (n in 1:50) {dens = sum(replicate(10000,{pMsw = pM; pMs = sample(1:50,n); pMsw[pMs] = pMsw[pMs[order(-(1:n))]]; mm1 = sample((1:50)[-c(4,15,43,11,18)],1); matches = (pMsw[c(4,15,43,mm1)] == pM[c(4,15,43,mm1)]); prod(matches==c(FALSE,TRUE,TRUE,FALSE))})); probdist_dens_new = rbind(probdist_dens_new,data.frame(n=n,dens=dens)) }
plot(probdist_dens_new,type="l")
with(probdist_dens_new,sum(n*dens)/sum(dens))
#[1] 26.03425
The mean seems to have shifted to a higher value, in agreement with one more wrong sample.
But I confess that at this point I'm not really sure what I'm doing... this is more complicated than it seemed.
