Models of Evolution - Revision history

Gorm at 14:50, 11 March 2026

2026-03-11T14:50:06Z

← Older revision		Revision as of 16:50, 11 March 2026
Line 386:		Line 386:
	ggplot(dflong, aes(x=div_time, y=value, col=name)) +		ggplot(dflong, aes(x=div_time, y=value, col=name)) +
	geom_line() +		geom_line() +
			geom_point() +
	labs(x="Time since divergence (MY)",		labs(x="Time since divergence (MY)",
	y="Genetic distance (substitutions/site)")		y="Genetic distance (substitutions/site)")

Gorm: /* The Kimura 2 parameter model */

2026-03-11T13:48:40Z

The Kimura 2 parameter model

← Older revision		Revision as of 15:48, 11 March 2026
Line 198:		Line 198:
	You should have obtained values of approximately A = -1.909 and B = -0.1818. You can now plot the curves showing how the expected amount of transitions and transversions change as a function of the branch length (the actual amount of change):		You should have obtained values of approximately A = -1.909 and B = -0.1818. You can now plot the curves showing how the expected amount of transitions and transversions change as a function of the branch length (the actual amount of change):
	df = tibble(		df = tibble(
	t = seq(0, 40, 0.1),		t = seq(0, 30, 0.1),
	Transitions = 0.25-0.5exp(At)+0.25exp(Bt),		Transitions = 0.25-0.5exp(At)+0.25exp(Bt),
	Transversions = 0.5 - 0.5 * exp(B*t),		Transversions = 0.5 - 0.5 * exp(B*t),

Gorm at 13:45, 11 March 2026

2026-03-11T13:45:06Z

← Older revision		Revision as of 15:45, 11 March 2026
Line 172:		Line 172:
	The branch length t measured in this way is related to ordinary calendar time as follows: if <math>\tau</math> denotes calendar time and <math>\mu=a+2b</math> is the overall substitution rate measured in substitutions per site per unit time, then <math>t=\mu \tau</math>.		The branch length t measured in this way is related to ordinary calendar time as follows: if <math>\tau</math> denotes calendar time and <math>\mu=a+2b</math> is the overall substitution rate measured in substitutions per site per unit time, then <math>t=\mu \tau</math>.

	Thus branch length ~~is simply~~ elapsed time multiplied by the substitution rate.		Thus branch length t can be thought of as time measured on an evolutionary scale: elapsed time multiplied by the substitution rate.

	This scaling is useful because sequence data usually inform us about the '''product''' of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time. To estimate absolute times and rates separately, the analysis must include calibration information, for instance known sampling dates for viral sequences at the tips.)		This scaling is useful because sequence data usually inform us about the '''product''' of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years.

	The formulas above therefore describe the expected observable ~~number~~ of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.		Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time. To estimate absolute times and rates separately, the analysis must include calibration information, for instance known sampling dates for viral sequences at the tips.

			The formulas above therefore describe the expected observable fractions of sites showing transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.

	For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of branch length <math>(t=1)</math> corresponds to <math>10^{9}</math> years.		For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of branch length <math>(t=1)</math> corresponds to <math>10^{9}</math> years.

Gorm at 13:36, 11 March 2026

2026-03-11T13:36:32Z

← Older revision		Revision as of 15:36, 11 March 2026
Line 174:		Line 174:
	Thus branch length is simply elapsed time multiplied by the substitution rate.		Thus branch length is simply elapsed time multiplied by the substitution rate.

	This scaling is useful because sequence data usually inform us about the '''product''' of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time.)		This scaling is useful because sequence data usually inform us about the '''product''' of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time. To estimate absolute times and rates separately, the analysis must include calibration information, for instance known sampling dates for viral sequences at the tips.)

	The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.		The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.

Gorm at 13:32, 11 March 2026

2026-03-11T13:32:43Z

← Older revision		Revision as of 15:32, 11 March 2026
Line 174:		Line 174:
	Thus branch length is simply elapsed time multiplied by the substitution rate.		Thus branch length is simply elapsed time multiplied by the substitution rate.

	This scaling is useful because sequence data usually inform us about the product of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time.)		This scaling is useful because sequence data usually inform us about the '''product''' of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time.)

	The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.		The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.

Gorm at 11:57, 11 March 2026

2026-03-11T11:57:15Z

← Older revision		Revision as of 13:57, 11 March 2026
Line 174:		Line 174:
	Thus branch length is simply elapsed time multiplied by the substitution rate.		Thus branch length is simply elapsed time multiplied by the substitution rate.

	This scaling is useful because sequence data usually inform us about the product of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (~~The exception is~~ when we fit ~~a concrete~~ clock~~-model to the data, where rate matrices are expressed~~ in ~~substitutions per site per time~~, ~~and~~ branch lengths in calendar time~~. We will do this using the software BEAST later~~ in ~~this course~~.)		This scaling is useful because sequence data usually inform us about the product of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (Later in the course, when we fit clock models in BEAST, branch lengths will instead be measured in calendar time and rates in substitutions per site per unit time.)

	The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.		The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.
Line 180:		Line 180:
	For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of branch length <math>(t=1)</math> corresponds to <math>10^{9}</math> years.		For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of branch length <math>(t=1)</math> corresponds to <math>10^{9}</math> years.

	In these scaled units, ~~branch length and time are numerically identical~~: a branch of length t corresponds to an expected t ~~substitutions~~ per site.		In these scaled units, the numerical value of t is the expected number of substitutions per site: a branch of length t corresponds to an expected t substitution events per site.

	----		----

Gorm at 11:55, 11 March 2026

2026-03-11T11:55:32Z

← Older revision		Revision as of 13:55, 11 March 2026
Line 148:		Line 148:
	We will now use RStudio to explore some features of evolution occurring according to this model.		We will now use RStudio to explore some features of evolution occurring according to this model.

	It can be shown that, under the K2P model, the ~~chance~~ of ~~observing~~ a transition and a transversion respectively ~~depends~~ on R and the branch length t in the following way:		It can be shown that, under the K2P model, the expected fractions of sites showing a transition and a transversion respectively depend on R and the branch length t in the following way:

	<math>		<math>
Line 168:		Line 168:
	</math>		</math>

	In the equations above, t is the branch length ~~expressed not in units of time but~~ as the ~~genetic distance measured in~~ expected ~~mutations per site. (Note: mutations mean how many actual events have happened, and this is different from substitutions, which is the~~ number of observable changes). The expected genetic distance t is directly proportional to the time elapsed. Specifically: if <math>\tau</math> denotes calendar time in some unit and <math>\mu=a+2b</math> is the overall substitution ~~rate (substitutions~~ per site ~~per unit calendar time,~~ measured in ~~the same unit), then <math>t=\mu \tau</math> (note how time factors out of this equation, and t now measures expected mutations per site, rather than~~ years or generations~~). Writing the equations in terms of t is convenient because the overall substitution rate is then 1 in these scaled units~~.		In the equations above, t is the branch length, measured as the expected number of substitution events per site. It is therefore not measured directly in years or generations.

	For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of t (1 ~~mutation per site~~) corresponds to ~~a clock-based branch length of~~ <math>~~\mu=~~10^{9}</math> years.		The branch length t measured in this way is related to ordinary calendar time as follows: if <math>\tau</math> denotes calendar time and <math>\mu=a+2b</math> is the overall substitution rate measured in substitutions per site per unit time, then <math>t=\mu \tau</math>.

			Thus branch length is simply elapsed time multiplied by the substitution rate.

			This scaling is useful because sequence data usually inform us about the product of rate and time, not about rate and time separately. A branch can become long either because much time has passed or because evolution proceeds quickly. For this reason, phylogenetic branch lengths are often expressed in substitutions per site rather than in years. (The exception is when we fit a concrete clock-model to the data, where rate matrices are expressed in substitutions per site per time, and branch lengths in calendar time. We will do this using the software BEAST later in this course.)

			The formulas above therefore describe the expected observable number of transitions and transversions as a function of branch length t, i.e. as a function of the total amount of evolutionary change that has accumulated.

			For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of branch length <math>(t=1)</math> corresponds to <math>10^{9}</math> years.

	In these scaled units, branch length and time are numerically identical: a branch of length t corresponds to an expected t substitutions per site.		In these scaled units, branch length and time are numerically identical: a branch of length t corresponds to an expected t substitutions per site.
Line 225:		Line 233:

	'''Question: '''Include the plots in your response. Which statement best describes the effect of increasing R?		'''Question: '''Include the plots in your response. Which statement best describes the effect of increasing R?
	:# It makes the transition curve rise faster and ~~peak earlier relative to the transversion curve~~.		:# It makes the transition curve rise faster and reach its maximum at a shorter branch length.
	:# It makes transitions dominate for longer at short distances, but also saturate more strongly.		:# It makes transitions dominate for longer at short distances, but also saturate more strongly.
	:# It has no effect on the asymptotic total observed difference.		:# It has no effect on the asymptotic total observed difference.
Line 381:		Line 389:
	: We have here plotted the total difference, the observed transitional and transversional difference, as well as the JC- and K2P-corrected distances as a function of divergence times.		: We have here plotted the total difference, the observed transitional and transversional difference, as well as the JC- and K2P-corrected distances as a function of divergence times.

	'''Question: '''Compare ~~the~~ p-distance, JC distance, and K2P distance as divergence time ~~increases~~. Why do ~~the~~ corrected distances increasingly exceed ~~the observed~~ p-distance? ~~Compare the~~ JC and ~~K2P corrections. What~~ does ~~the difference suggest~~ about ~~the relative rates of transitions and transversions~~ in this mitochondrial data set?		'''Question: '''Compare how p-distance, JC distance, and K2P distance change with divergence time. Why do both corrected distances increasingly exceed p-distance? Why is K2P consistently larger than JC, and what does this imply about transition/transversion bias in this mitochondrial data set?

	----		----
Line 399:		Line 407:
	: This shows the observed ratio between transitions and transversions as a function of divergence times.		: This shows the observed ratio between transitions and transversions as a function of divergence times.

	'''Question: '''Describe how the observed transition/transversion ratio changes with divergence time (and compare to the theoretical ~~behaviour you saw in~~ Question 6). Why is this ~~pattern~~ expected under the K2P model~~, and what~~ does it suggest about the true underlying transition/transversion rate ratio R in these mitochondrial sequences?		'''Question: '''Describe how the observed transition/transversion ratio P/Q changes with divergence time, and compare this with the theoretical behavior from Question 6. Why is this decline expected under the K2P model? What does it suggest about the true underlying transition/transversion rate ratio R in these mitochondrial sequences?

Gorm at 11:49, 11 March 2026

2026-03-11T11:49:10Z

← Older revision		Revision as of 13:49, 11 March 2026
Line 212:		Line 212:
	: Third, note that the expected amount of transitional differences first rise rapidly and then decline slowly to an equilibrium value of 0.25. Transversional differences rise slowly to an equilibrium value of 0.5. The equilibrium values are determined by the fact that when sufficient time has passed sequence similarities will essentially be random; since there are twice as many possible transversions as transitions, these will in the end make up two thirds of all observed changes. Early on, before this stage is reached, the much higher rate of transitions will cause them to make up the vast majority of all observed changes, and only after considerable time has elapsed will the transversions catch up.		: Third, note that the expected amount of transitional differences first rise rapidly and then decline slowly to an equilibrium value of 0.25. Transversional differences rise slowly to an equilibrium value of 0.5. The equilibrium values are determined by the fact that when sufficient time has passed sequence similarities will essentially be random; since there are twice as many possible transversions as transitions, these will in the end make up two thirds of all observed changes. Early on, before this stage is reached, the much higher rate of transitions will cause them to make up the vast majority of all observed changes, and only after considerable time has elapsed will the transversions catch up.

	'''Question: '''Estimate from the plot at approximately what branch length the observed ~~number of transitions~~ reaches its maximum. ~~Explain why~~ the ~~number of~~ observed ~~transitions~~ then ~~decreases~~, even though ~~transitions still~~ occur ~~in the evolutionary process.~~		'''Question: '''Estimate from the plot at approximately what branch length the observed transition fraction reaches its maximum. Why does the observed transition fraction then decline, even though transition events continue to occur?

	----		----

Gorm at 11:07, 11 March 2026

2026-03-11T11:07:38Z

← Older revision		Revision as of 13:07, 11 March 2026
Line 4:		Line 4:
	== Overview ==		== Overview ==

	In this exercise we will explore ~~a number of~~ different, but closely related, models of evolution. Using such models it is possible to estimate the number of unseen mutational events and thereby obtain genetic distances that have been corrected for superimposed substitutions. It is, however, important to realize that these corrections are based on the assumption that we observe approximately the expected amount of change - if, for instance, 20 mutational events end up leading to no observable changes then it is impossible to guess the actual amount of change regardless of which correctional scheme we employ. Using more and longer sequences helps ensuring that the observed change is closer to the expected change, so the correction is more likely to be accurate with adequate amounts of data. The same models also play an important role in phylogenetic reconstruction based on maximum likelihood and Bayesian techniques.		In this exercise we will explore two different, but closely related, models of evolution. Using such models it is possible to estimate the number of unseen mutational events and thereby obtain genetic distances that have been corrected for superimposed substitutions. It is, however, important to realize that these corrections are based on the assumption that we observe approximately the expected amount of change - if, for instance, 20 mutational events end up leading to no observable changes then it is impossible to guess the actual amount of change regardless of which correctional scheme we employ. Using more and longer sequences helps ensuring that the observed change is closer to the expected change, so the correction is more likely to be accurate with adequate amounts of data. The same models also play an important role in phylogenetic reconstruction based on maximum likelihood and Bayesian techniques.

	== Getting started ==		== Getting started ==
Line 98:		Line 98:
	:<math>d_{JC}=-\frac{3}{4}ln\left( 1 - \frac{4}{3}p \right)</math>		:<math>d_{JC}=-\frac{3}{4}ln\left( 1 - \frac{4}{3}p \right)</math>

	Note that this correction will only work if the observed difference is approximately as expected. Consider this: In the dice-rolling simulation we found that if there has been 0.67 changes per site then the expected observed difference is 0.44, because some sites will experience more than one substitution. However, as you saw in the simulation, the actual observed difference can, for random reasons, be different from the expected 0.44 (say, 0.33 or 0.58). As an extreme example: if all mutations act on a single nucleotide site, then only 1 mutation is visible, regardless of the actual number of mutations. More generally: If the observed difference is not the same as the expected observed difference, then we will obviously also get the wrong estimate of the real distance after correcting for multiple substitutions.		Note that this correction will only work if the observed difference is approximately as expected. Consider this: In the dice-rolling simulation we found that if there has been 0.67 changes per site then the expected observed difference is 0.44, because some sites will experience more than one substitution (superimposed substitutions, this phenomenon is known as "saturation"). However, as you saw in the simulation, the actual observed difference can, for random reasons, be different from the expected 0.44 (say, 0.33 or 0.58). As an extreme example: if all mutations act on a single nucleotide site, then only 1 mutation is visible, regardless of the actual number of mutations. More generally: If the observed difference is not the same as the expected observed difference, then we will obviously also get the wrong estimate of the real distance after correcting for multiple substitutions.

	We will now plot (estimated) real change as a function of observed difference (this is the inverse of what you did before). In RStudio enter:		We will now plot (estimated) real change as a function of observed difference (this is the inverse of what you did before). In RStudio enter:
Line 148:		Line 148:
	We will now use RStudio to explore some features of evolution occurring according to this model.		We will now use RStudio to explore some features of evolution occurring according to this model.

	It can be shown that, under the K2P model, the chance of observing a transition and a transversion respectively depends on R and t in the following way:		It can be shown that, under the K2P model, the chance of observing a transition and a transversion respectively depends on R and the branch length t in the following way:

	<math>		<math>
Line 168:		Line 168:
	</math>		</math>

	In the equations above, time ~~has been rescaled by~~ the ~~overall substitution rate~~. If <math>\tau</math> denotes calendar time and <math>\mu=a+2b</math> is the overall substitution rate (substitutions per site per unit calendar time), then ~~we define~~ <math>t=\mu \tau</math>.		In the equations above, t is the branch length expressed not in units of time but as the genetic distance measured in expected mutations per site. (Note: mutations mean how many actual events have happened, and this is different from substitutions, which is the number of observable changes). The expected genetic distance t is directly proportional to the time elapsed. Specifically: if <math>\tau</math> denotes calendar time in some unit and <math>\mu=a+2b</math> is the overall substitution rate (substitutions per site per unit calendar time, measured in the same unit), then <math>t=\mu \tau</math> (note how time factors out of this equation, and t now measures expected mutations per site, rather than years or generations). Writing the equations in terms of t is convenient because the overall substitution rate is then 1 in these scaled units.

	Thus t measures expected substitutions per site, rather than years or generations. In other words, t is the branch length under this model. Writing the equations in terms of t is convenient because the overall substitution rate is then 1 in these scaled units.		For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of t (1 mutation per site) corresponds to a clock-based branch length of <math>\mu=10^{9}</math> years.

	For example: if <math>\mu=10^{-9}</math> substitutions per site per year, then 1 unit of t corresponds to <math>\mu=10^{9}</math> years.

	In these scaled units, branch length and time are numerically identical: a branch of length t corresponds to an expected t substitutions per site.		In these scaled units, branch length and time are numerically identical: a branch of length t corresponds to an expected t substitutions per site.
Line 362:		Line 360:

	'''Plot distances'''		'''Plot distances'''

	In RStudio enter:		In RStudio enter:
	df = read_table("titv.data", col_names = FALSE, skip=1)		df = read_table("titv.data", col_names = FALSE, skip=1)
Line 387:		Line 386:
	'''Question 13'''		'''Question 13'''

	'''Compute and plot empirical transition/transversion ratio'''		'''Compute and plot empirical transition/transversion rate ratio'''

	In RStudio enter:		In RStudio enter:
	df = df %>% mutate(titv = transitions/transversions)		df = df %>% mutate(titv = transitions/transversions)

Gorm at 10:12, 11 March 2026

2026-03-11T10:12:11Z

← Older revision		Revision as of 12:12, 11 March 2026
Line 399:		Line 399:
	: This shows the observed ratio between transitions and transversions as a function of divergence times.		: This shows the observed ratio between transitions and transversions as a function of divergence times.

	'''Question: '''Describe how the observed transition/transversion ratio changes with divergence time. Why is this pattern expected under the K2P model, and what does it suggest about the true underlying transition/transversion rate ratio R in these mitochondrial sequences?		'''Question: '''Describe how the observed transition/transversion ratio changes with divergence time (and compare to the theoretical behaviour you saw in Question 6). Why is this pattern expected under the K2P model, and what does it suggest about the true underlying transition/transversion rate ratio R in these mitochondrial sequences?