Everything you could want to know

- Infections first observed in 1970s: most post-transmission hepatitis were “non-A, non-B hepatitis” (Harvey Alter, NIH)
- Virus first isolated in 1988, team led by Michael Houghton
- Alter recipient of 2000 Albert Lasker Award for Clinical Medical Research for his work leading to the discovery of the virus causing hepatitis C infection (with Michael Houghton)
- transmissible via fluids (sex and IV drug use common modes)
- About 80% of cases lead to chronic (long-term) infection
**Estimated 150-200 million are chronically infected worldwide! That’s roughly 3% of the world’s population.**- Infects
*hepatocytes*(liver cells); long-term effects include cirrhosis, hepatocellular carcinoma (liver cancer). - No vaccine
- First therapies mid-1990s: interferon, interferon augmented with antiviral ribavirin
- 2011: introduction of highly effective antivirals (protease inhibitor therapies, different mechanisms of action)
- Current standard of treatment: combination of antivirals
**Latest**: antiviral sofosbuvir (trade name Sovaldi) in combination with one of a few other antivirals (e.g. declatasvir) = cure in a high proportion of cases! But in the US, Sovaldi costs $1000 a pill…

**HCV life cycle** (Nature Reviews Gastroenterology & Hepatology 6, 403-411 (July 2009))

*Data from* Neumann et al. (1990).

**Patient 1e (left)**

[ * Data : hide , shown as table , shown as CSV shown as QR * ]

```
# time (days), viral load (log10 per mL)
#
# HCV viral load data following interferon therapy
# Patient 1E from Neumann et al. (1998)
# http://www.jstor.org/stable/2897780
#
0.0000,100
0.0833,92
0.1666,140
0.2916,110
0.4166,98
0.5833,41
0.7916,12
1.0000,7
1.2083,7.5
1.4167,7
2,8
3,2.4
4,3.9
5,2.6
7,1.2
9,0.4
11,0.45
14,0.05
```

**Patient 3d (right)**

[ * Data : hide , shown as table , shown as CSV shown as QR * ]

```
# time (days), viral load (log10 per mL)
#
# HCV viral load data following interferon therapy
# Patient 3D from Neumann et al. (1998)
# http://www.jstor.org/stable/2897780
#
0,38
0.08333,44
0.29167,45
0.41667,6
0.58333,3.8
0.79167,1.7
1,1
1.2083,0.65
1.4167,0.5
2,0.15
3,0.2
5,0.05
7,0.05
9,0.01
```

Let \(I\) represent the number of infected hepatocyte cells at time \(t\) and \(V\) the number of virus, **in a single mL**. Which can be written in ODE form as \[\begin{align}
\frac{dI}{dt} &= \beta(T) V - \delta I
\\
\frac{dV}{dt} &= p I - c V
\end{align}\]

Having a model is great, but just the start of the story. In order to actually use this model, we have to get estimates for the various parameters we have introduced, \(c\), \(p\), \(\delta\), and \(\beta(T)\). That’s allot harder to do than you might first expect when dealing with a micro-organism inside a human body. In fact, before the discovery of some effective treatments, nobody knew how to do it.

At this point, we will apply a very powerful-but-sometimes-incorrect idea to simplify this system of equations – the separation of time-scales. The idea, in rough terms, is that different parts of large systems often change at different rates. Some parts change quickly, while other parts change slowly. If these speed differences are large enough, then we can sometimes use them to simplify the dynamics of a system. When a component of a system oscillates very quickly, we can try to replace that component with it’s average behavior. And when a component of a system changes very slowly, we can try to replace that component with its steady-state behavior. These are core concepts in the field of applied mathematics known as “asymptotic approximation theory”. While relatively simple to grasp, the mathematics of applying these ideas can be quite complicated, and sometimes leads to serious mistakes if carelessly applied.

Here in our HCV model, though, the idea turns out to be straight-forward and very useful. In long-term, chronic infection, prior to the start of treatment, we expect the number of target cells, infected cells, and virus particles to be roughly constant, day-to-day. Once treatment starts, the number of infected cells and virus will change rapidly, but the number of target cells will not change very quickly because the processes in the human body that control target-cell creation are relatively slow processes. Thus, we can hypothesize that the number of target cells will be at approximately the same steady-state value for the duration of our experiments. At steady-state prior to treatment, the time-derivatives of both the number of infected cells and the amount of virus should vanish, so \[\begin{align} 0 &= \beta(T) V - \delta I, \\ 0 &= p I - c V. \end{align}\] Solving these equations, we find \(\beta(T) = c\delta/p\). Substituting back in, our potentially non-linear system for viral dynamics becomes a linear 2x2 system. \[\begin{align*} \frac{dI}{dt}&= \frac{c\delta}{p} V - \delta I, \\ \frac{dV}{dt}&= pI-cV. \end{align*}\]

Estimate parameters using data following initiation of therapy. But… what does interferon therapy do?

- Decrease infection rate (\(c \delta/ p\))?
- Interferon with viral production (p)?

As it turns out, interferon therapy blocks the ability of cells to produce virus particles – the more effective it becomes, the fewer virus particles get produced. We can parameterize interferon treatment efficacy, then, as the fraction of viral production blocked, \(\epsilon\) with the natural constraint that \(0\leq \epsilon\leq 1\). Introducing this into our linear system, then, \[\begin{align} \frac{dI}{dt} &= \frac{c\delta}{p} V - \delta I \\ \frac{dV}{dt} &= (1-\epsilon)pI - c V \end{align}\] or in matrix form, \[\begin{gather} \frac{d}{dt}\left(\begin{array}{c}I\\V\end{array}\right) = \left(\begin{array}{cc}-\delta & c\delta/p \\(1-\epsilon)p&-c\end{array}\right) \left(\begin{array}{c}I\\V\end{array}\right) \end{gather}\]

We can solve this system to get \(V(t) = C _ 1 e^{-\lambda _ 1 t} + C _ 2 e^{-\lambda _ 2 t}\), where \(\lambda_1\), \(\lambda_2\) are the eigenvalues of the matrix \(\left(\begin{array}{cc}-\delta & c\delta/p \\(1-\epsilon)p&-c\end{array}\right)\),

\[\lambda_{+,-}=-\frac{c+\delta\pm\sqrt{(c-\delta)^2+4c\delta(1-\epsilon)}}{2}\]

An effective therapy should have \(\epsilon\) near one, making \(1-\epsilon\) small, in which case \[\lambda_{+,-}=-c,\;\;-\delta.\]

Note that our data appears **biphasic**: two phases of exponential decay (y-axis is on a log scale). Our modeling predicts that the slopes of the decay phases will give rough estimates for the infected cell death rate \(\delta\) and the viral clearance rate \(c\):

[Show code]

```
Pat1e=loadtxt('NeumannEtAl1998_Pat1e.csv',delimiter=',')
# Start by showing data for viral decay phase
figure(1)
plot(Pat1e[:,0], log10(10**4*Pat1e[:,1]), 'b.', markersize=15)
axis([0, 1.1, 4.5, 6.5])
xlabel('Time (Days)',fontsize=24)
ylabel('log$_{10}$[Viral load (per mL)]',fontsize=24)
legend(['Data'],loc=1,fontsize=20)
xticks(fontsize=20)
yticks(fontsize=20)
plt.tight_layout()
savefig('FirstPhaseData_1e.png')
# do linear fit for viral clearance
coeff = polyfit(Pat1e[4:8,0], log(10**4*Pat1e[4:8,1]), 1)
t=linspace(0,14,1000)
Phase1=exp(coeff[1])*exp(coeff[0]*t)
plot(Pat1e[4:8,0], log10(10**4*Pat1e[4:8,1]), 'g.', markersize=15)
plot(t,log10(Phase1),'g--',linewidth=2)
legend(['Data','First phase data', 'First phase fit'],loc=1,fontsize=20)
savefig('FirstPhaseData_EstimateC_1e.png')
# Now show second phase data
axis([0, 14.2, 2, 6.5])
savefig('SecondPhaseData_1e.png')
# do linear fit for infected cell death rate
coeff2 = polyfit(Pat1e[8:,0], log(10**4*Pat1e[8:,1]), 1)
t=linspace(0,14,1000)
Phase1=exp(coeff2[1])*exp(coeff2[0]*t)
plot(Pat1e[8:,0], log10(10**4*Pat1e[8:,1]), 'r.', markersize=15)
plot(t,log10(Phase1),'r--',linewidth=2)
legend(['Data','First phase data', 'First phase fit','Second phase data', 'Second phase fit'],loc=1,fontsize=20)
savefig('SecondPhaseData_EstimateDelta_1e.png')
### Plot lambda1, lambda2 as a function of epsilon assuming c=0.39 and delta=2.57 as a function of epsilon
figure(2)
c=-coeff2[0]
delta=-coeff[0]
epsilon=linspace(0,1,100)
lambda1=-(c+delta)/2+sqrt((c+delta)**2-4*c*delta*epsilon)/2
lambda2=-(c+delta)/2-sqrt((c+delta)**2-4*c*delta*epsilon)/2
plot(epsilon,lambda1,'r-',linewidth=2)
#plot(epsilon,-0.39*ones_like(epsilon),'r--',linewidth=2)
plot(epsilon,lambda2,'b-',linewidth=2)
#plot(epsilon,-2.57*ones_like(epsilon),'b--',linewidth=2)
xlabel('IFN efficacy $\epsilon$',fontsize=24)
ylabel('Eigenvalues',fontsize=24)
legend(['$\lambda_1$; $\lambda_1(1) = -\delta$','$\lambda_2$; $\lambda_2(1) = -c$'],loc=2,fontsize=20)
xticks(fontsize=20)
yticks(fontsize=20)
plt.tight_layout()
savefig('Eigenvalues_1e.png')
```

For patient 1e the estimates are \(\delta=0.36\) per day and \(c=4.57\) per day. Based on our model we estimate that in patient 1e, infected cells die after an average of \(1/\delta\approx 3\) days.

**Model prediction: how many virus per day?** We can use our estimates to predict how many virus are produced in an infected individual per day. In our model that quantity is \(pI\). In the presence of therapy, \(\dot{V}=0=pI-cV\), so the quantity we’re looking for is \(=cV\). We prefer \(cV\) because we have an observations of \(V\): the average viral load before therapy. For patient 1e, the average is \(V_0=5.5\times 10^5\) virus per mL. Thus we predict that \(2.6\times 10^9\) virus produced per day per mL. Assuming 15L extracellular fluid (average human) that’s **\(\mathbf{3.9\times 10^{10}}\) HCV virions produced per day in an infected individual!**