Testing whether small deletions are under-represented within fibres containing multiple deletion species

mtDNASpecies.R

# Testing whether small deletions are under-represented within fibres containing multiple deletion species

# Assume 200 mtDNA molecules per fibre section
Nassume = 200

# P7, P15 and P16 are the proportions of smaller mtDNA species in fibres with two or more mtDNA species
P7 = c(0.6,0.2,0.49,0.33,0.57,0.75,0.47,0.29,0.27,0.23,0.51,0.54)
N7 = rep(Nassume,length(P7))
Ndel7 = c(2,2,3,2,2,2,2,3,2,2,2,2)
P15 = c(0.39,0.73,0.37,0.17,0.43,0.53,0.54,0.57)
N15 = rep(Nassume,length(P15))
Ndel15 = c(3,2,2,2,2,2,2,2)
P16 = c(0.23, 0.69, 0.86, 0.28, 0.25, 0.42, 0.67, 0.8, 0.79, 0.86)
N16 = rep(Nassume,length(P16))
Ndel16 = c(2,2,2,2,2,2,2,2,2,2)

getP = function(count,N,ndel = 2){
  return(binom.test(count,N,p=1/ndel,alternative="less")$p.value)
}

# Table 3 only includes proportions, not absolute numbers, so assume number of mtDNA molecules observed 
# in each fibre (N) is 200 for the moment
p7 = mapply(getP,round(P7*N7),N=N7,ndel=Ndel7)
p15 = mapply(getP,round(P15*N15),N=N15,ndel=Ndel15)
p16 = mapply(getP,round(P16*N16),N=N16,ndel=Ndel16)

# Correct p-values for multiple testing
pvals = p.adjust(c(p7,p15,p16),method="fdr")

# Count number of times smallest deletion (biggest molecule) makes up less than 50% of the mtDNA population
sum(pvals<0.05)
length(pvals)