# simple simulation of recurrent gene mutations

Sorry if this question is too basic for the forum. I am trying to determine the number of recurrently mutated genes expected by chance given an equal likelihood of mutation for all genes. There are about 20,000 genes and I have ~400 mutations spread across ~100 individuals. I want to determine the probability of the same gene being mutated in multiple individuals. I have tried the following simulation. I'm surprised that I'm always getting the same number of recurrences. Am I making some mistake?

  var = 400
individuals = 100
genes = 2e4

mutated_genes <- sample(rep(c(1,0), c(var,genes*individuals-var)))

sample_mat <- matrix(mutated_genes, ncol = individuals, nrow = genes)

trials <- 1e4

recurrences <- vector(mode = "list", length = trials)

for (i in seq_len(trials)){
# print(i)

new_rows <- sample(nrow(sample_mat))
new_cols <- sample(ncol(sample_mat))
new_mat <- sample_mat[new_rows, new_cols]
recurrences[[i]] <- rowSums(new_mat)
}

n_recurrences <- sapply(recurrences, function(x) sum(x > 1))

summary(n_recurrences)
$$$$
`
• Please say more about what you mean by "I have ~400 mutations spread across ~100 individuals." Is that 400 mutations total, for about 4 mutations on average per individual? Or is it that each individual has 400 mutations? In the latter case, is that each individual having exactly 400, or with some variability among individuals? I'm having some trouble trying to parse your code. It might help if you could annotate each line of the code with what you intend that line to accomplish.
– EdM
Feb 21 '21 at 22:58

If what you have is $$100$$ individuals each having $$400$$ mutations in one of $$20,000$$ genes, with all genes having the same probability of mutation, then for a quick analysis you don't need a simulation.

Overall, there are $$400 \times100=40,000$$ mutations. With $$20,000$$ genes, that comes out to 2 mutations per gene on the average over all $$400$$ individuals. If that's the number you were always coming up with, your simulations were making sense. There are lots of fun things you can play with in terms of the distributions of mutations among individuals and among genes in your set of $$400$$ samples, but the basic point is that in this scenario you should be surprised if you don't find many genes mutated in more than 1 sample.

In practice, as I'm sure you know, things aren't that simple. For the benefit of others reading this, here are some of the issues:

If you're starting with tumor samples, you are more likely to find cancer driver genes mutated because the tissues from which you are sampling have already been selected for being cancerous. For example the TP53 gene, the "guardian of the genome," is very frequently mutated in tumors.

Even beyond genes specifically driving tumor development, different genes have different probabilities of being mutated. That depends on gene size, gene transcriptional activity, and the timing of gene replication during the cell cycle. Lawrence et al. analyzed those situations in detail, working out a way to correct for those influences to determine whether particular genes are mutated at a higher rate than chance in a set of tumor samples.