Haddox
commented
2 years ago Here is a simplified version of the above idea, including a new term that accounts for functional constraint. It is inspired by equations 1-3 from a recent phydms paper, where now, $x$ and $y$ represent codons. When a mutation event happens in antigen, we randomly select the mutation it corresponds to by drawing from the following probability distribution:
$$P{r,xy} = \frac{Q{xy} \times F_{r,xy}}{Z}$$
where:
- $x$ is the wildtype codon
- $r$ is the codon number
- $y$ is the mutant codon (single nucleotide change away)
- $Q_{r,xy}$ is the rate of the mutation according to the HKY85 model
- $F_{r,xy}$ is the rate that the mutation fixes, as determined by amino-acid preferences
- $Z$ is a partition function that ensures that all probabilities sum to one
Specifically, we define $F{r, xy}$ in the same way as the phydms paper:
$$F{r,xy} = \frac{\ln{[ (\pi{r,y}/\pi{r,x}})^\beta ]}{1-(\pi{r,y}/\pi{r,x})^\beta}$$
where:
- $\pi_{r,x}$ is the preference for the amino acid encoded by codon $x$ at site $r$
- $\pi_{r,y}$ is the preference for the amino acid encoded by codon $y$ at site $r$
- $\beta$ is the stringency parameter
One way that the above idea violates the HYK model is that viruses can only get a maximum of one mutation per gene per day, since there is a maximum of one mutation event per day. Though, the mutation rate in antigen is quite low ( $10^{-4}$ mutations per host per day), so this might not be very much of a problem.
Erick pointed out that a more accurate approach would be to take a gene, cycle over each nucleotide, and apply the rate matrix (e.g., $Q{r,xy} \times F{r,xy}$) to each site. This would solve the above violation, as well as make things more interpretable. Though we don't know how much this would slow down the simulation since we'd need a for-loop over each site in the sequence.
matsen
In
antigen, there is a mutation rate $\mu$. In the original version ofantigen, which did not model viral sequences, mutations changed a virus's antigenic phenotype. In our updated version, which does model sequences, each mutation also results in a single nucleotide change.When a mutation event occurs in
antigen, we need to decide how to mutate the sequence. Our current strategy involves two basic steps. First, we randomly choose a site to mutate in the gene (e.g., C12). Second, we randomly choose a nucleotide-level mutation given the wildtype nucleotide and a pre-specified transition/transversion ratio $\kappa$. This approach is similar to nucleotide-mutation models used in phylogenetics, but differs in a few ways.The goal of this issue is to update the mutation model in
antigenso that mutations occur at rates proportional to ones from the HKY85 mutation model. In this model, the rate at which nucleotide $x$ mutates to nucleotide $y$ is:$$Q_{xy} = \begin{cases} \kappa\pi_y & \text{if } x \text{ and } y \text{ differ by a transition}\ \pi_y & \text{if } x \text{ and } y \text{ differ by a transversion}\ \end{cases}$$
where $\pi_y$ gives the expected frequency of nucleotide $y$ in the absence of selection. Our current mutation model in
antigendiffers from the HKY85 model in two ways. First, it does not account for $\pi_y$ parameters in choosing $y$. Second, by choosing a random site to mutate in the first step from above, it ignores that a site's mutation rate depends on the wildtype nucleotide $x$.Below is an idea for how to update the mutation model in
antigen. A simple idea is to replace the single $\mu$ parameter with $Q_{xy}$ parameters. However, inantigen, $\mu$ accounts for both mutation and selection. The virus population in a host is represented by a single virus with a defined antigenic phenotype and gene sequence. According to the literature, mutations seem to rarely fix in hosts. Thus, $\mu$ could be interpreted as the rate at which a mutation first occurs, then increases to an appreciable frequency, and then is involved in a transmission event. Since it is difficult to model each component of this process, the below idea retains the concept of $\mu$, but ensures that mutations occur at rates proportional to the HKY85 model. Specifically, when a mutation event occurs, the new model uses the following workflow to determine how to mutate the nucleotide sequence.Step 1: choose a site to mutate: In the HKY85 model, the probability that a site mutates depends on its wildtype nucleotide. Specifically, the rate at which a nucleotide $x$ mutates to any of the other three nucleotides is:
$$R_x = \sumy{Q{xy}}$$
where the expression sums over each possible mutant nucleotide. Under this model, if there is a single nucleotide mutation somewhere in a gene sequence, then the probability that the mutation involves a wildtype nucleotide of $x$ is:
$$P_x = \frac{R_x n_x}{\sum_i{R_i n_i}}$$
where $n_x$ is the number of counts of nucleotide $x$ in the gene sequence, and the expression in the denominator sums over each of the four nucleotides.
We will mimic this process in the following way. When there is a mutation event in
antigen, we will first randomly select one of the four nucleotides to be the wildtype nucleotide in the mutation, with selection probabilities of $P_x$. Then, we will make a list of all sites in the gene sequence with the randomly selected wildtype nucleotide (e.g., C12, C24, C28, etc.), and randomly choose one of those sites to mutate (e.g., C12). This list could be stored in a dictionary and updated for computational efficiency.Step 2: choose a mutant nucleotide: Once we have chosen a site to mutate, we will randomly select a mutant nucleotide $y$ with probabilities proportional to mutation rates $Q_{xy}$.
Using these two steps, $Q{xy}$ values in
antigenshould be proportional to $Q{xy}$ values in HKY85, unless I am mistaken. What do you all think? Any other ideas?