Skip to content

Instantly share code, notes, and snippets.

@sathomas
Last active August 29, 2015 14:16
Show Gist options
  • Star 0 You must be signed in to star a gist
  • Fork 0 You must be signed in to fork a gist
  • Save sathomas/a0c05c989f5466d786f3 to your computer and use it in GitHub Desktop.
Save sathomas/a0c05c989f5466d786f3 to your computer and use it in GitHub Desktop.
Mathematical Models of HIV, Part 2

In part 1 we looked at a basic mathematical model for HIV within a human host. This visualization continues that discussion by looking at one factor that limits the initial proliferation of the virus.

In part 1 we looked at a basic mathematical model for HIV within a human host. That model consists of differential equations using three variables:

  • x : the concentration of healthy CD4+ T cells (a favorite target of the virus),
  • y : the concentration of CD4+ T cells infected by the virus, and
  • v : the concentration of free virus particles circulating in the host.

The resulting model is

$ [ [ dot x ] , [ dot y ], [ dot v ] ] = [ [ lambda - d x - beta x v ], [ beta x v - a y ], [ k y - u v ] ] $

Even though its equations are very simple, there are no analytical solutions for the three variables as functions of time. We cannot solve for exact equations, so we turn to other techniques for insight from the model.

Equilibrium Points

Equilibrium points tells us where the system variables become constant. If the system reaches an equilibrium point, it remains in that state indefinitely. Part 1 notes that our model has two equilibrium points, and it discusses the disease-free equilibrium where both y and v are zero (and x = λ/d.) Under some conditions, that equilibrium point is stable—if the system gets close enough, then internal dynamics take over and drive it the rest of the way to the equilibrium. Unfortunately, those conditions do not appear to exist for HIV.

Endemic Equilibrium

In the case of HIV it is the second equilibrium point, the endemic equilibrium point, that is stable. The following graph reproduces the chart from part 1, but with a single curve representing today’s best estimates for parameters corresponding to HIV. After several weeks, the amount of virus in the patient’s blood settles at a relatively constant level.

1101001k10k100k1M10M100MFree Virus Particles020406080100120140Days Since Infection

Target Cell Limited

Another obvious quality of the chart above is the extremely high peak in the first few days after infection. The virus load reaches almost 1 billion particles before falling back down. The three equations of our model do not (yet) account for the body’s immune system, so the decline cannot be due to the actions of the body’s defenses. It is, instead, an intrinsic property of the model.

To see the reason for the decline, we can consider not only the amount of virus, but also the amount of healthy and infected CD4+ T cells in the host. In our model, it is the healthy cells that are the virus’s target. If there are no target cells available, then the virus cannot increase its numbers. In fact, that is exactly what the model predicts. The headline graph in this article focusses on the first three weeks of infection. It shows the number of free virus particles as a red line, and it also shows the percentage of the body’s CD4+ T cells that remain healthy. As the graph shows, at about the third day the percentage of healthy cells drops precipitously, almost to zero. At exactly that point the virus runs out of targets, and its numbers begin their decline. In the initial phase of an HIV infection, the virus’s own efficiency is the main factor limiting its growth.

In the next installment we’ll use a less traditional approach to visualize the target limiting effect.

References

Disclaimer

I’m neither a virologist nor a mathematician, just a geek that likes math. Undoubtedly some of what I’ve said above is wrong. When you notice it, please let me know, and I’ll do my best to fix it.

<!DOCTYPE html>
<html>
<head>
<meta charset='utf-8'>
<title>Mathematical Models of HIV, Part 2</title>
<link href='http://fonts.googleapis.com/css?family=Varela'
rel='stylesheet' type='text/css'>
<style>
body { font-family: Varela,sans-serif; color: #444; }
</style>
</head>
<body>
<script src='http://d3js.org/d3.v3.min.js'></script>
<script>
// Parameters for the simulated system. Explanations
// for these parameters may be found in the associated
// text.
var d = 0.01, // death rate of healthy cells
a = 1.00, // death rate of infected cells
u = 23.0, // death rate of free virus particles
k = 1000, // rate of virus generation
lambda = 100000, // rate of healthy cell generation
beta = 2e-8; // viral efficiency
// The derived parameters that the simulation actually
// used.
var alpha = a / d, // α = a/d
epsilon = d / u, // ε = d/u
R0 = (k * lambda * beta) / (a * d * u); // R0 = kλβ/dau
// Convenient ratios to convert from system state to
// meaningful physical quantities.
var ul2person = 5e6, // adults have about 5L of blood
t2days = 1 / d; // convert system time to days
// The system we're simulating is three differential
// equations. To keep the function definitions clean,
// we get the system parameter values from globals
// instead of function parameters.
function dx_dt(x,y,v) { return 1 - x - R0 * x * v; };
function dy_dt(x,y,v) { return R0 * x * v - alpha * y; };
function dv_dt(x,y,v) { return ( alpha * y - v ) / epsilon; };
// The time increment used in the simation.
var dt = 0.0001;
// Here's the function that actually simulates the system.
// It's inputs are the initial value for v, the number of
// points to generate, and the number of time steps between
// each point.
var simulate = function(v0, numPoints, numSteps) {
var result = [];
// The simulation begins with 100% healthy
// cells (and no infected cells.
var prev = { x: 1, y:0, v: v0, t: 0 };
for (var i=0; i<numPoints; i++) {
// At the beginning of each t-cycle, add a new
// data point to the results. Convert from the
// derived units used in the simulation to
// meaningful physical quantities when we do so.
result.push({
t: prev.t * t2days,
v: prev.v * ul2person,
x: prev.x * lambda / d * ul2person,
y: prev.y * lambda / d * ul2person,
r: prev.x / (prev.x + prev.y)
});
// Now quickly step through all the individual
// delta-t values for the current t-cycle.
for (var j=0; j<numSteps; j++) {
// Make sure none of the values go negative
// (and that v stays above the minimum).
var x = prev.x + dx_dt(prev.x, prev.y, prev.v) * dt,
y = prev.y + dy_dt(prev.x, prev.y, prev.v) * dt,
v = prev.v + dv_dt(prev.x, prev.y, prev.v) * dt;
// Make sure none of the values go negative.
prev.x = Math.max(0,x);
prev.y = Math.max(0,y);
prev.v = Math.max(1/ul2person,v);
prev.t += dt;
}
}
return result;
};
// Define the dimensions of the visualization. The
// width and height dimensons are conventional for
// visualizations displayed on http://jsDataV.is.
var margin = {top: 30, right: 60, bottom: 48, left: 60},
width = 636 - margin.left - margin.right,
height = 476 - margin.top - margin.bottom;
// Create the SVG stage for the visualization and
// define its dimensions. Web browsers don't require
// the full SVG attributes, but adding them makes it
// easier to extract the SVG from a web page and
// insert it into a stand-alone file.
var svg = d3.select("body").append("svg")
.attr("width", width + margin.left + margin.right)
.attr("height", height + margin.top + margin.bottom)
.attr("viewBox", "0 0 " +
(width + margin.left + margin.right) + " " +
(height + margin.top + margin.bottom) )
.attr("version", "1.1")
.attr("xmlns", "http://www.w3.org/2000/svg");
// Within the SVG container, add a group element (<g>)
// that can be transformed via a translation to account
// for the margins. This element is restricted to the
// clipping path.
var g = svg.append("g")
.attr("transform", "translate(" + margin.left +
"," + margin.top + ")");
// Define the scales that map a data value to a
// position on the vertical axis and to a time
// on the horizontal axis.
var v = d3.scale.log(),
t = d3.scale.linear(),
r = d3.scale.linear();
// Define a convenience function to create a line on
// the chart. The horizontal axis (which D3 refers
// to as `x`) are the time values, and the vertical
// axis (which D3 refers to as `y`) are the x-values.
// The result is that `line` is function that, when
// passed a selection with an associated array of
// data points, returns an SVG path whose coordinates
// match the t- and x-scales of the chart.
var line = d3.svg.line()
.x(function(d) { return t(d.t); })
.y(function(d) { return v(d.v); });
// Do the same to create an area for the fraction
// of healthy cells.
var area = d3.svg.area()
.x(function(d) { return t(d.t); })
.y0(height)
.y1(function(d) { return r(d.r); });
var datasets = [];
datasets.push(simulate(1/ul2person, width, 4));
// Set the scales
v.domain([1, 1e9]).range([height,0]);
r.domain([0,1]).range([height,0]);
t.domain([0, datasets[0][datasets[0].length-1].t])
.range([0, width]);
// Define a function that will create the axes
// when passed a data selection.
var tAxis = d3.svg.axis()
.scale(t)
.tickSize(5, 0)
.orient("bottom");
var vAxis = d3.svg.axis()
.scale(v)
.orient("right")
.ticks(5,"s");
var rAxis = d3.svg.axis()
.scale(r)
.tickFormat(d3.format(".0%"))
.ticks(5)
.orient("left");
// Create a selection for the data set
// and add the data to it.
var areas = g.selectAll(".area")
.data(datasets);
areas.enter().append("path")
.classed("area", true)
.attr("fill", "#7EBD00")
.attr("d", area);
var lines = g.selectAll(".line")
.data(datasets);
lines.enter().append("path")
.classed("line", true)
.attr("fill", "none")
.attr("stroke-width", "2px")
.attr("stroke", "#ca0000")
.attr("d", line);
// Add the t-axis.
svg.append("g")
.attr("class", "t axis")
.attr("transform", "translate("+margin.left+","+(height+margin.top)+")")
.call(tAxis)
.append("text")
.attr("x", width/2)
.attr("y", 36)
.attr("text-anchor", "middle")
.text("Days Since Infection");
// Add the r-axis.
var rgroup = svg.append("g")
.attr("class", "r axis")
.attr("transform", "translate("+margin.left+","+margin.top+")")
.call(rAxis);
rgroup.append("text")
.attr("x", -50)
.attr("y", -16)
.text("Target Cells (Healthy vs. Total)");
rgroup.append("rect")
.attr("x", 156)
.attr("y", -30)
.attr("width", 18)
.attr("height", 18)
.attr("fill", "#7EBD00");
// Add the v-axis.
var vgroup = svg.append("g")
.attr("class", "v axis")
.attr("transform", "translate("+(width+margin.left)+","+margin.top+")")
.call(vAxis);
vgroup.append("text")
.attr("text-anchor", "end")
.attr("x", margin.right)
.attr("y", -16)
.text("Free Virus Particles");
// Style the axis. We could do this with CSS, but
// by using JavaScript we ensure that the styles
// are embedded in the SVG itself. That's helpful
// if we ever want to extract the SVG from the web
// page as an image.
svg.selectAll(".axis line, .axis path")
.attr("fill", "none")
.attr("stroke", "#bbbbbb")
.attr("stroke-width", "2px")
.attr("shape-rendering", "crispEdges");
svg.selectAll(".axis text")
.attr("fill", "#444")
.attr("font-size", "14");
svg.selectAll(".axis .tick line")
.attr("stroke", "#d0d0d0")
.attr("stroke-width", "1");
// Append v legend after styling to avoid
// overwriting.
vgroup.append("line")
.attr("x1",-90)
.attr("y1",-20)
.attr("x2",-72)
.attr("y2",-20)
.attr("stroke","#ca0000")
.attr("stroke-width",2);
</script>
</body>
</html>
http://jsdatav.is/img/thumbnails/hiv2.png
Sign up for free to join this conversation on GitHub. Already have an account? Sign in to comment