An Intracellular Model of Hepatitis B Viral Infection: An In Silico Platform for Comparing Therapeutic Strategies

Hepatitis B virus (HBV) is a major focus of antiviral research worldwide. The International Coalition to Eliminate HBV, together with the World Health Organisation (WHO), have prioritised the search for a cure, with the goal of eliminating deaths from viral hepatitis by 2030. We present here a comprehensive model of intracellular HBV infection dynamics that includes all molecular processes currently targeted by drugs and agrees well with the observed outcomes of several clinical studies. The model reveals previously unsuspected kinetic behaviour in the formation of sub-viral particles, which could lead to a better understanding of the immune responses to infection. It also enables rapid comparative assessment of the impact of different treatment options and their potential synergies as combination therapies. A comparison of available and currently developed treatment options reveals that combinations of multiple capsid assembly inhibitors perform best.

and are located at the 3 and 5 ends of the minus strand, are removed. Then, a cellular 23 DNA polymerase, such as DNA polymerase κ, and either ligase I or ligase III attach and Step 2-Epigenetics of cccDNA: An important factor in liver disease caused by HBV 27 infections is the cccDNA. This molecule is organised into mini-chromosomes with host cell 28 histones and potentially other host and viral proteins [6,7]. Therefore, transcription from 29 cccDNA is subject to epigenetic regulation, with numerous options for dynamic epigenetic 30 control of cccDNA transcriptional activity [8,9,10]. In the absence of the X protein (HBx), Since the level of surface protein is higher than the HBV DNA level, we assume the 49 binding rates that lead to the production of the LS and S mRNAs to be n Ld times larger [5], 50 and b 1 = b 2 = b 3 = b. The parameters c 1 , c 2 , ...c 5 are elongation rates.

51
Step 4-Viral protein synthesis: All viral mRNAs have to be exported from the nucleus 52 to the cytoplasm in order to be translated [2]. Once a free ribosome successfully binds to 53 mRNA at rate d, a protein is produced and the ribosome releases the mRNA [15]. The 54 pgRNA is bicistronic and produces the C and P proteins. Thus, when the ribosome binds 55 to pgRNA, it can produce either the C or P protein. This has been modelled as follows: where e 1 and e 2 are translation rates. As production of a virion needs more C proteins than 57 polymerase, we assume that the pg:Rib 2 complex is produced at an n c -fold higher rate than 58 pg:Rib 1 . where e 3 , e 4 , ..., e 7 are translation rates. We assume the binding rates that lead to the pro-60 duction of surface proteins to be n Lm times larger [5]. Similarly to pgRNA, S mRNA is 61 also bicistronic, and we assume the S:Rib 2 complex to be produced at an n s -fold higher rate 62 because S protein is the most common protein in the viral envelopes and is known to result 63 in production of multiple non-viral particles [16]. 64 Note, in the protein production reactions, that these proteins are produced as monomers. 65 Subsequently, they bind to form dimers, which act as building blocks. However, for simplicity, 66 we assume in this model that these proteins are already in dimeric form, as removing this 67 step does not affect model behaviour.

68
The following equations describe the degradation of mRNAs and proteins: These degradation rates have been computed based on the half-lives of these mRNAs 69 and proteins. 70 Step 5-Capsid assembly initiation: Assembly of HBV depends on the formation of a 71 ribonucleoprotein (RNP) complex formed from P and pgRNA [2,17]. A short structured 72 RNA signal located at the 5 end of the pgRNA, called ε, was identified as the RNA packaging 73 signal mediating the packaging of pgRNA into the nucleocapsid (NC) [18,19]. ε was later 74 found to be recognised specifically by the P protein, not C, and P and pgRNA packaging 75 are mutually dependent [20,21]. The formation of this RNP complex depends on cellular 76 factors, which include the heat shock protein 90 (Hsp90), and also require ATP hydrolysis 77 and the function of p23, an identified chaperone partner for Hsp90 [3,17]. p23 can bind to 78 the polymerase independently of Hsp90. The binding of Hsp90 and p23 to the polymerase is 79 not sufficient to enable the polymerase to bind ε , but the Hsp90-p23 interaction is critical 80 for the formation of the viral RNP complex (binding of P to the pgRNA) [17]. Since Hsp90 81 and p23 are cellular proteins, whose production is independent of the virus, we assume that 82 they are available in sufficient quantities and are therefore not rate limiting. Thus, in the 83 mathematical model, we only represent the binding of pgRNA and P with rate g, with the 84 impact of these different proteins contained implicitly in the rate: Step 6a-Capsid assembly: The RNP complex is ready to be encapsidated. There are 86 Figure S1: All the possible reactions of the RNP complex (R) with C proteins.
Once the complex RC 1 C 2 C 3 is formed, more Cs are recruited to build the full capsid 94 according to the following reactions: where κ is an average binding rate of C to the growing capsid shell, and pgNC denotes a 96 fully formed capsid containing pgRNA, called the immature nucleocapsid.

97
Step 6b-Empty capsid assembly: It has been discovered recently that infected cells also 98 produce empty virions, i.e., empty viral capsids enclosed by surface proteins [14]. Creation 99 of an empty virion is made of two steps: first, the nucleation step where C proteins bind 100 together at a slower rate [23,24]. We model this step as follows: where nuc is the number of core proteins that are needed for creating the nucleation complex 102 and k n is the forward rate for the nucleation reactions.

103
The second step is the elongation step where the remaining C proteins assemble at the 104 faster rate κ 2 to form an empty nucleocapsid (emNC): In this model, we assume that empty capsids are made of 120 C proteins following a 106 T = 4 architecture, ignoring the T = 3 particles that occur in only 5% of cases [22].

107
Step 7-Reverse transcription: The reverse transcription of the pgRNA starts at this 108 stage. In ∼ 90% of cases, this process leads to the synthesis of rcDNA, and in the remaining 109 10%, to double-stranded linear DNA (dslDNA) [25]. Therefore, we have where rcNC and dslNC indicate nucleocapsids containing rcDNA and dslDNA, respectively.

111
These are called the mature nucleocapsids.

112
Step 8-Intracellular cccDNA amplification: rcNC and dslNC can deliver their contents  Thus, we have the following equations: where t is the export rate.
As the production rate of dslNC is much lower than that of rcNC (less than 10%) and 128 the integration occurs once in any ∼ 10 5 − 10 6 infected cells, and since the NHEJ is error 129 prone, many of cccDNAs produced from dslDNA are functionally defective, and we do not 130 integrate this step into the model and only consider the rcNCs.

131
Step 9a-Secretion of complete virions: When the level of the L protein is sufficiently Step 9b-Secretion of RNA-containing particles: It has been discovered recently that Step 9c-Secretion of empty virions: Like complete virions, we assume that the ratio of

162
Step 9d: Secretion of filaments: Tubular filaments are a form of SVP that are also called contain 120 surface proteins on average. Therefore, we have: where t 3 is the secretion rate of filamentous SVP. SVPs are released in much higher numbers than infectious virions, in order to simplify our 181 model, we assume that as soon as 48 S proteins are assembled, they are released as a spherical 182 particle, rather than creating larger filaments that are then divided into several particles.

183
Although this assumption is a simplification, it does not affect the overall outcome of the 184 model, as filaments keep budding spheres at a fast rate, which can be modelled by choosing 185 the right value for the budding rate.

186
It has been reported that a small proportion of L protein results in the production and where the impact of L protein is added to the initiation reaction, κ is the binding rate of S 191 proteins, and t 2 is the secretion rate of spherical SVPs. 193 In this section, we focus on the current and future treatment options for chronic hepatitis 194 B. Figure S14 shows all the treatment options and the process that they attack. Below, we 195 discuss details for each drug and introduce our mathematical models for their impacts. Thus, in the model, we consider that IFN-α binding silences the cccDNA and IFN-αincreases 204 the degradation rate of cccDNA. Hence, the additional reactions have the form:

S2 Effects of Drug Treatments
where IFN is the amount of IFN-α, ψ is the efficacy of this drug in degrading of the cccDNA, 206 and I bind and I unbind are the binding and unbinding rates of IFN-α, respectively. which cannot bind to the pgRNA. Therefore, we have the following reactions: where k bind and k unbind indicate the binding and unbinding rates of GA, respectively.  Since this drug (CAM) can bind to core proteins and alter their assembly pathway, we 254 add the following reactions regarding binding and unbinding of CAMs to core proteins: where C bind and C unbind are the binding and unbinding rates of this drug to C protein, 256 respectively.

257
The C:CAM complexes can bind together to form capsids. If the drug is class I, they 258 form morphologically intact capsids. Therefore, 120 of them bind together, at a faster rate 259 than C-C interactions (κ), to form empty capsids, and we assume that these capsids are where nuc is the number of core proteins that are needed for composing the nucleating 262 structure [23, 60], k nc is the forward rate for the nucleation reactions, and κ 1 is the forward 263 rate for the elongation reactions.

264
As suggested in [54], CAMs can bind to rcNCs and trigger their disassembly followed by 265 degradation of the rcDNA. We also assume that capsid proteins are partially assembled and 266 therefore not available for buildup of new capsids; hence, we only include a reaction that 267 results in CAM becoming available after dissociation, as follows: where k dis indicates the disassembly rate.

269
If the drug is class II, its binding to C causes the formation of aberrant capsids. In this 270 case, for simplicity, we assume that aberrant capsids contain 120 C:CAM on average, and 271 we model them in the same way, but we assume that they are degraded rather than released, 272 as the effect on the overall infection dynamics is equivalent in both cases.

273
Note that we focus the model on a more efficient version of CAMs, i.e., JNJ-827, which 274 is class I and has a smaller half-maximal effective concentration.
where siRNA is the number of small interfering RNAs, R bind is the binding rate of the 312 interfering RNA to the pgRNA, and k cleave is the cleaving rate, which leads to the degradation 313 of pgRNA.

314
Inhibitors of viral and SVP release: These drugs target the viral and SVP egress, i.e., 315 these drugs reduce the rates t, t 1 , t 2 , and t 3 . Let us assume that η be the efficacy of these is one-fourth of the product lifetime, and the activation and deactivation rates are equal. 334 We have used the same assumption. As the half-lives of proteins are considered to be one are typically found at 1000-100,000-fold and 100-fold higher levels than complete virions, 349 respectively, and RNA-containing particles are 100-1000-fold lower than complete virions.

350
Based on these facts, the other parameters are chosen in a way to see these dynamics (Table   351 S1), which has been shown in Figure 2c-f. For example, we assume that the budding rate of 352 RNA-containing particles is 1000 times smaller than that of complete virions (t 1 = t/1000).

353
Although the yield of empty virions is 100-fold higher than that of complete virions (note that 354 their assembly is faster because they do not need to encapsidate a viral genome), we assume 355 that they have the same budding rate as complete virions. As SVP yield is approximately 356 10 | | | | | 100-fold higher than that of empty virions, and since, similarly to empty virions, they do 357 not need to encapsidate a genome, we assume that their budding rate is 100 times higher 358 than that of empty virions (t 2 = 100t). Note that these are our approximations, which are 359 consistent with in vivo observations [14].

360
As shown in Figure S3, choosing the initial number of free ribosomes equal to 10,000 361 provides a more realistic intracellular dynamic, and this value is a reasonable choice based , k B is the Boltzmann constant, and T is the temperature in Kelvin [85].  It has been observed that increasing the binding rates of a drug does not have a significant 381 effect on the outcome of the treatment, perhaps because this will also increase the unbinding 382 rate, as the unbinding over binding rate is constant (dissociation constant). Thus, we assume 383 that the binding rate of IFN is equal to b (binding rate of RNA Pol II for pgRNA), and the 384 binding rate of NA to polymerase is equal to g (binding rate of polymerase and pgRNA).

385
The binding rate of siRNA molecules to the pgRNA is assumed to be equal to d (binding 386 rate of ribosome to pgRNA) with K d = 0.1 nM, and the cleaving rate k cleave = 10 hour −1

387
[86]. We also consider the PS-targeting drugs and CAMs to have the same binding rate as C 388 protein to PSs and k dis = 10 × ln(2)/24 = 0.29, which means that the half-life of rcNC:CAM 389 is 2.4 hours ( Figure S8b). In addition, we assume that the binding rate of C:CAMs is f 390 times faster than the binding rate of C proteins for the production of an empty virion, i.e.,

391
κ nc = f κ n and κ 1 = f κ. We assume that f = 10. However, changing 10 to a smaller or 392 bigger value does not have a significant effect on the outcome of CAMs ( Figure S8a). For 393 the inhibitor of HBsAg release, we assume that its efficacy is equal to 95%.  Sens: Sensitivity analysis was performed. * These values are approximated based on Figure 7A in [72]. * * Measured in an animal system. * * * Sensitivity analysis was performed, and this value was chosen to get SVPs at a 10,000-fold higher level than complete virions.      Table S1 and α = 1 hour −1 . Figure S7: The number of complete virions released with RNA interference treatment and parameter values from Table S1 as a function of the binding affinity (K d ) and the cleaving rate (k cleave ) of this drug. The concentration of the drug is 50 molecules.
16 Figure S8: The number of complete virions released with CAM treatment and parameter values from Table S1. In (a), k dis = 0.29 and concentration of CAM is 20 molecules. In (b), f = 10 and C bind = 24.     Table S1. The grey dashed line shows the start of treatment and the shaded area around each curve indicates the 95% CI. Figure S13: The HBV life cycle.
Step 0: Viral entry is mediated by the sodium-taurocholate co-transporting polypeptide (NTCP) receptor, which allows the nucleocapsid (rcNC) containing relaxed circular DNA (rcDNA) to be released into the cell via endocytosis. This step is not included explicitly in the model.
Step 1: After attachment of the rcNC to the nucleus, it delivers the rcDNA. The host DNA repair factors convert the rcDNA into covalently closed circular DNA (cccDNA).
Step 2: X protein, which gets produced by the cccDNA over the course of infection, promotes the de-silencing of the cccDNA and blocks its silencing.
Step 3: The cccDNA is used as a template by RNA Polymerase II (RNA Pol II) to synthesise viral RNAs, including the pgRNA; LS, S, PreC, and X mRNAs. The pgRNA encodes both the core (C) and polymerase (P) proteins; X mRNA the X protein; LS mRNA the L surface protein; S mRNA the M and S surface proteins; and PreC mRNA the precore (PreC) protein.
Step 4: Translation of these mRNAs by ribosomes (Rib) leading to synthesis of viral proteins.
Step 5: The pgRNA and P form a 1:1 ribonucleoprotein (RNP) complex, which is assembly competent.
Step 6a: Encapsidation of the RNP complex by C proteins to form a nucleocapsid containing pgRNA:P (pgNC), also termed the immature nucleocapsid.
Step 6b: Assembly of C proteins into an empty nucleocapsid (emNC).
Step 7: The reverse transcription of pgRNA by P within the pgNC, resulting in the conversion of pgNC into a rcDNA containing nucleocapsid (rcNC), which is also termed the mature nucleocapsid.
Step 8: Recycling of the rcDNA from mature nucleocapsids to form more cccDNAs.
Step 9a: Envelopment of mature nucleocapsids within a membrane layer containing the surface proteins L, M, and S, leading to the secretion of a complete virion.
Step 9b: Envelopment of an immature nucleocapsid resulting in secretion of an RNA-containing particle.
Step 9c: Envelopment of an empty nucleocapsid and secretion of an empty virion.
Step 9d: L, M, and S proteins form empty filaments and filamentous subviral particles (SVP), via conventional tubular budding, into late endosomes and exit the cell.
Step 9e: S proteins assemble into octahedral spheres (spherical SVPs) and exit the cell.