Theoretical Explanation for the Rarity of Antibody-Dependent Enhancement of Infection (ADE) in COVID-19

Global vaccination against the SARS-CoV-2 virus has proved to be highly effective. However, the possibility of antibody-dependent enhancement of infection (ADE) upon vaccination remains underinvestigated. Here, we aimed to theoretically determine conditions for the occurrence of ADE in COVID-19. We developed a series of mathematical models of antibody response: model Ab—a model of antibody formation; model Cv—a model of infection spread in the body; and a complete model, which combines the two others. The models describe experimental data on SARS-CoV and SARS-CoV-2 infections in humans and cell cultures, including viral load dynamics, seroconversion times and antibody concentration kinetics. The modelling revealed that a significant proportion of macrophages can become infected only if they bind antibodies with high probability. Thus, a high probability of macrophage infection and a sufficient amount of pre-existing antibodies are necessary for the development of ADE in SARS-CoV-2 infection. However, from the point of view of the dynamics of pneumocyte infection, the two cases where the body has a high concentration of preexisting antibodies and a high probability of macrophage infection and where there is a low concentration of antibodies in the body and no macrophage infection are indistinguishable. This conclusion could explain the lack of confirmed ADE cases for COVID-19.


Introduction
Since 2019, coronavirus disease (COVID- 19) has spread all over the world. Severe forms of COVID-19 are associated with acute respiratory distress syndrome, pneumonia, renal failure and death [1,2]. Global vaccination is a path to the reduction of the disease's severity, as well as its spreading rate [3]. Nevertheless, among the general concerns for the development and utilization of vaccines is the possibility of antibody-dependent enhancement of infection (ADE) [4]. Besides vaccination, ADE is frequently caused by viruses that have high antigenic diversity and the ability to replicate in immune cells [5], such as coronaviruses, Ebola virus, dengue virus, human immunodeficiency viruses (HIV) and influenza [6][7][8][9][10]. ADE mostly occurs when a non-neutralizing or poorly neutralizing antibody (Ab) binds to a virus particle [11]. Normally, binding of the Ab to the Fc receptors (FcRs) on leukocytes leads to the destruction of the virus inside these cells [12]. However, in the case of ADE, Abs can release the virus and let it replicate in FcR-expressing cells [13]. Moreover, such infection of immune cells can cause adverse responses [11], such as enhanced cytokine release by the infected cells, which facilitate tissue damage [14].
ADE can originate from different features of (pre-existing) antibodies. For dengue infection, it was shown that both non-neutralizing antibodies and antibodies that are (a-c) model Ab: (a) approximation of experimental data for viral load dynamics [34]. The following parameter values were chosen (Equation (1)): ARNA = 798,292 RNA numbers, x0= 2.998 days, w1 = w2 = 2.069 days −1 , w3 = 2.38 days. (b) Changes in relative numbers of antibodies and dynamic of viral epitopes for an individual patient [34]. (c) Kinetics of naive B-cells, plasma cells and memory B-cells for the same patient. The following donor-specific parameters were used: δpl = 0.004 days −1 , πpl = 3•10 −5 (APS•days) −1 , βRNA = 264, w3 = 2.62 days, tnaive = 12.8 days, w = 1.04. (d-f) model Cv (d) A typical fit to the experimental data (black dots) on viral load in COVID-19 patients [1]. The day of seroconversion is marked with the grey dotted line and assumed to happen on day 6. The following donorspecific parameters were used: βinf =2. 21 The following average values of parameters for the given group of patients (6 patients with mild and 7 patients with severe course of disease) were obtained: βRNA = 532.38, (a-c) model Ab: (a) approximation of experimental data for viral load dynamics [34]. The following parameter values were chosen (Equation (1)): A RNA = 798,292 RNA numbers, x 0 = 2.998 days, w 1 = w 2 = 2.069 days −1 , w 3 = 2.38 days. (b) Changes in relative numbers of antibodies and dynamic of viral epitopes for an individual patient [34]. (c) Kinetics of naive B-cells, plasma cells and memory B-cells for the same patient. The following donor-specific parameters were used: δ pl = 0.004 days −1 , π pl = 3 × 10 −5 (APS·days) −1 , β RNA = 264, w 3 = 2.62 days, t naive = 12.8 days, w = 1.04. (d-f) model Cv (d) A typical fit to the experimental data (black dots) on viral load in COVID-19 patients [1]. The day of seroconversion is marked with the grey dotted line and assumed to happen on day 6. The following donor-specific parameters were used: β inf = 2.21 × 10 −6 (Virus·days) −1 , δ v = 255 days −1 , ac 0 = 2.17 × 10 −10 (cells·cells·days) −1 , Virus 0 = 1.32 × 10 7 particles; (e) The dependence of the viral response dynamics for different values of the probability on macrophage infection probability ρ. The following average values of parameters for the given group of patients (6 patients with mild and 7 patients with severe course of disease) were obtained: β RNA = 532.38, w 3 = 2.378 days, t naive = 10.5 days, δ pl = 0.166 days −1 , and π pl = 10 −4 (APS·days) −1 . Typical dynamics of naive B-cells, plasma cells and memory B-cells are given in Figure 1c, Figures S1 and S2. Altogether, model Ab correctly describes antibody response, B-cell differentiation and immune system-mediated virus elimination in response to the given viral load. However, this model does not take into account variabilities in the types of antibodies and the possible infection of macrophages.

The Infection Propagation Model (Model Cv) Suggests Possibilities of ADE Occurrence in COVID-19
The infection propagation model (model Cv) was developed to describe antibodydependent infection of macrophages. The model parameters were adjusted based on viral replication kinetics in HAE cells [35] ( Figure S4) and the viral load data from the sputum of seven COVID-19 patients from Germany [1] (Figures S3 and 1d). The model was able to describe viral load data for different patients upon variation of such parameters as pneumocyte infection rate β inf , the maximum virus-antibody complex formation constant α c0 , the rate of viral degradation δ v and the initial viral load Virus 0 . After approximation of seven experiments, the following average parameter values were obtained: β inf = 4.72 × 10 −6 (days·Virus) −1 , α c0 = 5.99 × 10 −4 days −1 , δ v = 207.6 days −1 and Virus 0 = 2.6 × 10 7 . Other parameters are given in Table S7. The dynamics of the viral load and the pneumocyte number for different ρ are given in Figure 1e,f.
We define ADE as a drastic drop in the concentration of healthy pneumocytes (more than 5% in comparison with the course of disease without macrophage infection). Model Cv could predict conditions for the occurrence of ADE for average antibody level. For this, we calculated the dependence of the viral load dynamics on the probability of macrophage infection (ρ) (Figure 1e). At ρ < 10 −3 , the concentration of viral particles was the same as in the absence of macrophage infection, and no ADE was observed, while for ρ > 10 −3 , the probability of pneumocyte survival decreased and the times of viral clearance increased, allowing the occurrence of ADE. However, we could not investigate the impact of antibody parameters on ADE occurrence using only model Cv because the existing experimental data describe antibodies in terms of optical density and we needed module model Ab to obtain the absolute concentrations of antibodies. Then, we estimated the impact of the concentration of pre-existing antibodies (Figure 2c,d). The model predicted that only the highest concentration (A 0 ≈ 10 13 particles) of antibodies could significantly change the course of disease; i.e., the number of surviving pneumocytes increased drastically (Figure 2d), while the viral load noticeably decreased ( Figure 2c). As discussed above, the dissociation constant (K d ) of the antigen-antibody complex could affect the course of the desease. To determine the extent of K d influence, we performed similar simulations for varied values of K d (Figure 2e-h). In the absence of pre-existing antibodies, Kd did not influence the viral load maximum, but larger Kd led to lower pneumocyte survival ( Figure 2f) and prolonged virus elimination (Figure 2e). In the presence of preexisting antibodies at the maximum concentration (A 0 ≈ 10 13 particles), the viral load decreased with time, all pneumocytes survived for small Kd values (Kd < 10 −9 M) and up to 20% of the pneumocytes died for large Kd values (Kd > 10 −9 M). Therefore, the existance of low-affinity pre-existing antibodies could lead to ADE. (c,d) Dynamics of free viral epitope (c) and pneumocytes (d) at different concentrations of pre-existing antibodies. Label Ab 0 denotes the concentration of pre-existing antibodies in lymphoid organs at the onset of disease. (e,f) Dynamics of free viral epitope (e) and pneumocytes (f) at various dissociation constants for the antigen-antibody complex for primary antibodies. (g,h) Dynamics of free viral epitope (g) and pneumocytes (h) at various dissociation constants for the antigen-antibody complex for pre-existing antibodies. The initial concentration of pre-existing antibodies equaled 10 13 particles (Ab lo = 10 13 ).

Theoretical Conditions for ADE Occurrence in COVID-19
The next question of interest was whether ADE could occur in primary infection.
To answer this question, we calculated the number of pneumocytes, the viral load and the number of life macrophages in the absence of pre-existing antibodies ( Figure 3). Our results showed that ADE could occur in the absence of pre-existing antibodies if the probability of macrophage infection was greater than 10 −3 . At ρ < 10 −3 , the concentration of viral particles was the same as in absence of macrophage infection. Upon an increase in the probability of macrophage infection, pneumocyte survival decreased and viral load clearance time increased.  (e,f) Dynamics of free viral epitope (e) and pneumocytes (f) at various dissociation constants for the antigen-antibody complex for primary antibodies. (g,h) Dynamics of free viral epitope (g) and pneumocytes (h) at various dissociation constants for the antigen-antibody complex for pre-existing antibodies. The initial concentration of pre-existing antibodies equaled 10 13 particles (Ablo = 10 13 ).

Theoretical Conditions for ADE Occurrence in COVID-19
The next question of interest was whether ADE could occur in primary infection. To answer this question, we calculated the number of pneumocytes, the viral load and the number of life macrophages in the absence of pre-existing antibodies ( Figure 3). Our results showed that ADE could occur in the absence of pre-existing antibodies if the probability of macrophage infection was greater than 10 −3 . At ρ < 10 −3 , the concentration of viral particles was the same as in absence of macrophage infection. Upon an increase in the probability of macrophage infection, pneumocyte survival decreased and viral load clearance time increased.    First of all, we investigated the influence of pre-existing antibody concentration on the course of disease (Figures 4 and S5) for a given probability of macropage infection (ρ). For ρ < 3 × 10 −4 , ADE was not observed (the percentage of dead pneumocytes did not exceed levels 3% higher than with the course of disease without macrophage infection; this situation occurred with extremely high antibody concentrations 10 14 < Ab old < 10 15 ). For ρ > 3 × 10 −4 , ADE occurred only if the concentration of pre-existing antibodis was high enough (Ab old > 10 10 ).
Secondly, we estimated the role of the antigen-pre-existing antibody complex dissociation constant (Kd) in possible ADE occurrence (Figures 5 and S6). Interestingly, the influence of Kd on the course of disease was depended on the macrophage infection probability ρ (see Texts S3 and S4 in the Supplementary Materials for details). For ρ < 3 × 10 −4 , the proportion of surviving pneumocytes increased upon a decrease in Kd. In this case, the maximal number of dead pneumocytes was no greater than their amount in the absence of antibodies and possible macrophage infection ( Figure 5). For ρ > 3 × 10 −4 , the highest level of pneumocyte death was observed at the lowest Kd. In all these cases, the degree of Kd impact on pneumocyte survival depended on the pre-existing antibody concentration ( Figure 5).

Discussion
The goal of our study was to theoretically determine the conditions for ADE occurrence in COVID-19.
First of all, we constructed an antibody response generation model that could describe antibody dynamics in COVID-19. We modified the mathematical model [36] for ADE in Dengue fever using parameters corresponding to SARS-CoV-2 infection to describe the observed viral load kinetics. The constructed model was able to describe the experimental viral load dynamics [1], seroconversion times and antibody concentrations [34]. We analyzed the influence of key parameters on the course of the disease (Figure 2). We found that an increase in the initial viral load induced an increase in the maximum viral load and caused an earlier rise and decline in the viral load. This result corresponds with the predictions of other SARS-CoV-2 mathematical models [37]. For other type of

Discussion
The goal of our study was to theoretically determine the conditions for ADE occurrence in COVID-19.
First of all, we constructed an antibody response generation model that could describe antibody dynamics in COVID-19. We modified the mathematical model [36] for ADE in Dengue fever using parameters corresponding to SARS-CoV-2 infection to describe the observed viral load kinetics. The constructed model was able to describe the experimental viral load dynamics [1], seroconversion times and antibody concentrations [34]. We analyzed the influence of key parameters on the course of the disease (Figure 2). We found that an increase in the initial viral load induced an increase in the maximum viral load and caused an earlier rise and decline in the viral load. This result corresponds with the predictions of other SARS-CoV-2 mathematical models [37]. For other type of virus (influenza), an increase in the initial viral load led to a time-shift of the peak without changing its magnitude [38]. This difference could be explained by differences in the assumed host cell behavior; specifically, in our model, there was no mechanism for pneumocyte recovery after infection, unlike other mathematical models for respiratory infections [38,39].
In this study, we assumed that ADE occurs if the amount of dead pneumocytes is higher in the presence of pre-existing antibodies. It should be noted that, in contrast to previously published models [36], in the mathematical model constructed ere the pneumocyte damage varied from patient to patient (10-90%) and never reached 100%.
Using the complete model, we found that the probability of macrophage infection regulated the impact of antibody concentration over the course of disease (Figures 4 and 5). At any probability of macrophage infection, the concentration of pre-existent antibodies played a key role in antibody-mediated enhancement of disease. This result corresponded with a previously predicted high impact of antibody concentration on ADE occurance [39], although in that work no explicit immune cells were included. Our calculations demonstrated that, with moderate macropage infection probabilities, ADE was possible only for the highest antibody concentrations (Figures 4 and S5). Similar results were obtained experimentally for dengue disease [40]. The lower the concentration of antibodies in the body was, the higher the probability of macrophage infection required in order to cause significant changes in the kinetics of the viral load and survival of pneumocytes. When the concentration of antibodies was high, the proportion of infected macrophages increased. However, in such cases, antibodies began to have a neutralizing effect. As a result, the number of dead pneumocytes in cases with possible infections of macrophages and a high antibody titer was no greater than the number of dead pneumocytes in cases where the level of antibodies was low and there was no antibody-mediated infection of macrophages. Thus, even if an antibody-dependent macrophage infection is detected, it is most likely that the lungs will not be more affected than in the course of the disease in the absence of ADE.
The antigen-antibody affinity can affect the initiation of ADE, but its influence on the course of the disease is controlled by the level of antibodies in the body. A marked decline in pneumocyte survival will only be observable when the binding affinity between the viral epitope and the antibody is very high. This suggestion is indirectly supported by a clinical case in which an acute myeloid leukemia patient might have suffered ADE shortly after being treated with high-affinity anti-COVID-19 therapy [41].
However, a situation in which the antigen-antibody affinity for the newly synthesised antibodies is lower than for pre-existing antibodies seems extremely unlikely, significantly lowering the probability of ADE after vaccination. Indeed, consistent with [42,43], it can now be assumed that anti-SARS-CoV-2 vaccination does not aggravate the disease course in COVID-19.
There is a set of experimental findings that were not taken into account by the developed model and which could have possibly made it underestimate the effect that non-neutralizing antibodies have on the course of disease. First, it is noteworthy that lung epithelial cells also express high levels of FcγRIIa, and non-neutralizing antibodies might facilitate their infection [44]. In addition, COVID-19 patients have been reported to have a strong IgG antibody response to the nucleocapsid protein, which cannot neutralize SARS-CoV-2 at all, resulting in delays in virus clearance and an increased severity of infection [44].
Another point to be discussed is the development of ADE upon treatment with convalescent plasma. Although the mathematical model proposed here cannot be directly applied in this case, some data can be obtained from such studies. Altogether, the model results are in accordance with the data from [45], where it was shown that, in a cell culture, donor and patient plasma could induce ADE, but its probability did not differ between plasma samples from convalescent, mild, moderate and severe patients. In contrast to this data, a study on the transfer of neutralization activity with COVID-19 convalescent transfusion [46] demonstrated the absence of ADE upon transfusion. Based on our calculations, we speculate that convalescent plasma transfer is not likely to change the course of disease, but other plasma transfusion risks, such as hemolytic transfusion reactions, anaphylactic reactions, transfusion-related acute lung injury and allergic reactions, might occur [47].

General Principles for the Construction of the Mathematical Models
The models contained several compartments, with the assumption of well-mixing of species inside each compartment (homogeneous models). The models were based on ordinary differential equations derived from the laws of chemical kinetics (the mass action law, Michaelis-Menten kinetics or Hill equation). The schemes of the models are given in Figure 1. Model parameters were taken from the literature data or estimated based on the available experimental data, as explained in the Supporting Information.
In order to describe the ADE phenomenon through mathematical modeling, we defined ADE in the text as a drastic drop in the concentration of healthy pneumocytes (more than 5% in comparison to the course of disease without macrophage infection).

Model Ab-The Antibody Response
Model Ab was constructed to describe interactions between the pathogen (a virus particle) and the antibody-production system (Figure 6a). It included induction of Ab synthesis in lymphoid organs, binding of the Abs to viral particles and Ab-dependent degradation of viruses. Model Ab comprised two different compartments: lungs and lymphoid organs. B-cell activation, proliferations, differentiation and antibody-production occur in lymphoid organs, whereas viral interaction with antigen-presenting cells and antibodies takes place in the lungs. We assumed that the virus-antibody complex formed between a single viral protein (epitope) and an antibody independently of steric interference between the proteins on the viral surface.
In model Ab, it was assumed that the number viral particles (Virus) was proportional to the number of SARS-CoV-2 RNA in patient samples. The kinetics of the SARS-CoV-2 RNA number in saliva was set analytically as: where F RNA is the concentration of SARS-CoV-2 RNA, t is the model time and A RNA , x 0 , w 1 , w 2 and w 3 are parameters required to better fit the experimental data [34]. Numerical values of these parameters are represented in Table S1. It was assumed that the parameter w 3 was specific for individual patients. Individual parameter values are represented in Table S6.
In the frame of model Ab, we assumed that each SARS-CoV-2 virus particle had exactly 25 S proteins on its surface and each protein on the virus could bind antibodies independently [48]. Thus, the number of viral particles was proportional to the sum of the viral epitope bound with the antibody (Va) and the free viral epitope (V f ). Subsequently, the number of free viral epitope could be described with the following formula: where β RNA is the constant of proportionality between the SARS-CoV-2 RNA number and the amount of free viral epitopes. It was assumed that this parameter implied the specificity of the viral kinetics in individual patients. Individual parameter values are represented in Table S2. tion: where Ablung is the concentration of antibodies in the lungs; kv1 and kv−1 are the epitopeantibody association and dissociation rates, respectively; cv1 is the maximal rate of clearance of viral protein-antibody complexes; and cv2 is the half-decay constant. In model Ab, it was assumed that the clearance rate of the epitope-antibody complex was approximated by the Michaelis-Menten formula, as described in the following equation: where Ab lung is the concentration of antibodies in the lungs; k v1 and k v−1 are the epitopeantibody association and dissociation rates, respectively; c v1 is the maximal rate of clearance of viral protein-antibody complexes; and c v2 is the half-decay constant. The next equations for the model describe the Ab production in the response to a viral load. Formation of the long-term immune response is a complex process governed by multiple factors, such as the number of antigen-presenting cells, T-cell help, the area of contact between T-and B-cells on the border of lymph nodes, etc. [48]. In order to avoid excessive complexity, we unified all of these factors through the introduction of antigen-presenting sites (APSs). APSs are represented for the lungs and lymphoid organs. Lung APSs can exist in two states: neutral (Ap) and activated (APS lung ). Activation of APSs occurred in the lungs when viral epitopes interacted with antigen-presenting cells. B-cells could be activated in the lymphoid organs. Transition of activated antigen-presenting sites in the lung (APS lung ) into activated antigen-presenting sites in lymphoid organs (APS lo ) occurred in accordance with the law of mass action. In the absence of a pathogen, all sites were assumed to be in a neutral state, Ap, and they could not affect B-cell functioning. The dynamics of neutral APSs and activated APSs in lungs and lymphoid organs can be described by the following equations: dAps lo dt = k s1 * Aps lung − k s−1 * Aps l0 − δ apm Aps lo (6) where χ b is the generation and degradation of neutral APSs in the absence of the virus (homeostasis rate of neutral APSs); Ap 0 and Ap are the initial and transient concentrations of neutral APSs; δ ap and δ apm are the rates of APS lung and APS lo deactivation, correspondingly; the terms β apm * Ap * V f /(ϕ + V f ) indicate the rate of APS lung activation, β apm is the maximal activation rate; and ϕ is the half-saturation constant for APS lung activation. ϕ was estimated from an existing model [49,50]. The terms k s1 and k s−1 denote rates of APS transition between two compartments: lungs and lymphoid organs. ν is the conversion coefficient reflecting the lymphoid organ:lung volume ratio. Activation of B-cells was assumed to occur in a separate compartment in lymphoid organs, such as the spleen or lymph nodes [48]. After interaction with activated sites (APS lo ), naive B-cells (B) proliferated and differentiated into short-lived IgM-producing plasma cells, long-lived IgG-producing plasma cells (Pl) and memory B-cells (Bm) (Figure 5a). We focused on B-cells, which produce high-affinity class-switching IgG-antibodies. The possibility of B-cell differentiation into short-lived IgM-producing plasma cells was also included in the model. In the frame of the model, we assumed that B-cell proliferation and differentiation into short-lived plasma cells occurred after their contact with APSs, whereas differentiation into long-lived plasma cells and memory B-cells occurred after a time delay. The complete equation for the concentration of B-cells is the following: Aps lo B − Aps lo B τ n * π pl + τ n * π bm + π ps (7) where B 0 and B are the initial and transient concentrations of B-cells; χ b is the B-cell generation and degradation in the absence of the virus (homeostasis rate of B-cells); π b is the naive B-cell proliferation rate (reflecting the increase in B-cell proliferation upon contact with APSs) [51]; the π parameters reflect probabilities of different B-cell fates upon contact with the APSs: π pl -into plasma cells, π bm -into memory cells and π ps -into short-lived Ig-M-producing plasma cells; τ n is used to introduce a time delay between the appearance of APSs and the initiation of naive B-cell differentiation into long-lived plasma cells and memory B-cells [48]; t 0 is the moment when APS lo becomes positive; and t naive is the lag-time required for B-cell differentiation. The next equation for model Ab described the dynamics of IgG-producing plasma cells (Pl): dPl dt = τ n π pl Aps lo B + τ mem π ml Aps lo B m − δ pl Pl (8) where δ pl is the rate of plasma cell degradation and τ mem reflects the time lag (t memory ) between the appearance of APSs (at t 0 ) and the initiation of memory B-cells differentiation and proliferation [48].
As well as the formation from B-cells and production of plasma cells, the memory B-cells were assumed to be capable of reproduction with saturation [48]. Thus, the memory cell dynamics were described by the following equation.
where k bm1 and k bm2 represent the B m growth rate and the lymphoid organ maximum capacity, respectively. The concentration of antibodies in the lymphoid organs and body fluids (A b ) could be described with the following equation: where α pl denotes the rate of IgG generation by long-living plasma cells [52] and δ a is the natural decay rate of antibodies (A lo and A lung ) [38,51]. The terms k a1 and k a−1 denote rates of antibody transition between two compartments, the lungs and lymphoid organs.
To compare the model predictions with experimental data on serum antibody levels [34], we introduced into the model Ab parameters ω and Ig0, which indicate the background level and a scaling coefficient that depended on the type of antibody test. Thus, the relative IgG equaled ln(Ab lo /Ig0+ω).

Model Cv-Propagation of Infection
Model Cv was constructed to describe the elimination of viral particles by alveolar macrophages, viral entry and replication in pneumocytes and antibody-mediated viral entry into macrophages (Figure 6b). Model Cv described events localized in a lung alveolus; therefore, the parameters differed from model Ab.
The antibodies were assumed to be produced by plasma cells, the concentration of which could be described by a simplified version of Equation (8): where κ b and µ B are the plasma cell production and degradation rates and α B is the rate of plasma cell proliferation upon contact with viral particles [36]. As antibodies are produced by plasma cells, in model Cv we assumed the concentration of antibodies to be proportional to the concentration of plasma cells and, in the following equations, Pl mod is used to indicate Ab concentration and Pl mod * Virus is used to indicate the concentration of virus-antibody complexes.
In model Cv, macrophages were assumed to be in two states: healthy (M A ) and infected (M*). Healthy macrophages engulfed viral particles without changing their concentration or engulfed virus-antibody complexes and became infected. In the model, it was assumed that the infected macrophages could recover [36] and both healthy and infected macrophages could degrade. With the described assumption, the equations for M A and M* were the following: where χ M is the macrophage generation and degradation in the absence of a virus (homeostasis rate of macrophages); δ I is the death rate of infected macrophages (δ I ≥ χ M ); σ is the rate of infected macrophage recovery; ρ denotes the probability of macrophages becoming infected; and α c is the rate for antibody-virus complex formation and subsequent engulfment, expressed as a step function: where α c0 is the maximal constant rate, t is the model time and t ser is the day of seroconversion. The function of the rate for antibody-virus complex formation and the corresponding parameters were chosen to allow approximation of antibody kinetics with model Ab. We used t ser * 2 to acquire a better approximation of the antibody kinetics. Healthy pneumocytes (P n ) were assumed to be the main virus-producing cell. They could engulf the virus, become infected (P i ) and, consequently, die [53]: where β inf is the rate of infection of pneumocytes and δ p is the death rate of infected pneumocytes. Thus, viral concentration was described with the following equation: where p is the rate of virus production in pneumocytes, α v is the rate of viral phagocytosis by the macrophages and δ v is the rate of macrophage-independent viral clearance. We assumed that the number of viral particles was proportional to the RNA concentration obtained by PCR in COVID-19 patients. Parameter αRNA was added to connect viral RNA concentration in sputum to the actual virus concentration in the lungs. It was individually chosen simultaneously for all patient data and selected to be in the range from 100 to 0.001.

Complete Model-A Combination of Model Ab and Model Cv
The complete model was designed to combine the antibody production from model Ab and the production of viruses and infection of macrophages from model Cv. The scheme of the model is given in Figure 6c. The required modifications of the models included scaling of the volumes and concentrations, addition of macrophage-independent degradation of both free and antibody-bound virus epitopes, conversion of viral particle concentrations (Viral) from model Cv into viral epitope concentration (V) from model Ab and direct inclusion of the virus-antibody complexes (Va) with the equations for M A and M*. Therefore, Equations (2), (3), (13) and (14) were modified as follows: Parameters π V * and ρ* were recalculated from p and ρ by dividing them by the number of S proteins on the viral particle. Here, we did not include viral production by macrophages in Equation (18). In some studies, SARS-CoV-2 replication in immune cells has been shown to be abortive, as no viable virions were produced [25][26][27].
In the complete model, the infection probability was linked to the features of the antibodies indirectly because it is currently assumed that a set of antigen-antibody interaction properties are suspected to cause ADE [15,54], and one of the most important properties of IgG antibodies is their capacity to engage and activate specific FcγR pathways, which does not depend on their affinity [16].

Model Integration, Analysis and Parameter Estimation
The systems of ordinary differential equations were integrated using the LSODA method [55,56] in COPASI software (http://www.copasi.org, accessed on 1 January 2020).
Estimations of model parameters were undertaken separately for model Ab and model Cv based on experimental data or previously published models, as described below and in the Supporting Information (Texts S1 and S2). Additional parameter adjustment for the complete model was based on experimental data [1,35,57].
For model Ab, values for naive B-cell proliferation, antibody secretion by plasma cells and antibody death rates were taken from published mathematical models and the literature [38,51,58], as indicated in Tables S1-S4. Parameters for the antigen-antibody complex formation were obtained from the average values for binding with monoclonal antibodies [59]. The ranges for time lags t memory and t naive were calculated from known data on lymphocyte physiology [60,61] and experimental results [34] (see Text S1 for the calculation details).
For model Cv, parameters were either fitted for personal patient data [1] or taken from [36] (Table S4).
For the complete model, the parameters π V * and δ p were estimated from the data on SARS-CoV-2 replication in human airway epithelial (HAE) cells [35]. The set of equations for HAE cell infection is represented in Text S3 in the Supplementary Materials. The parameters for the macrophage activity were calculated based on experimental data and published models [1,38,57], as described in Text S2 in the Supplementary Materials.
The analysis of the models consisted of time-course simulations, analysis of steady states performed using Newton's method [55] in COPASI and sensitivity analysis [55,62,63].
For the latter, we calculated sensitivities using the "one-at-a-time" approach and the following formula: where R i indicates the model response (Va, S a , A b , B, Lp, Bm at time 30 days for model Ab) and p j -s are the models parameters. The detailed results are represented in Table S5.