The LDL-HDL Profile Determines the Risk of Atherosclerosis: A Mathematical Model

Atherosclerosis, the leading death in the United State, is a disease in which a plaque builds up inside the arteries. As the plaque continues to grow, the shear force of the blood flow through the decreasing cross section of the lumen increases. This force may eventually cause rupture of the plaque, resulting in the formation of thrombus, and possibly heart attack. It has long been recognized that the formation of a plaque relates to the cholesterol concentration in the blood. For example, individuals with LDL above 190 mg/dL and HDL below 40 mg/dL are at high risk, while individuals with LDL below 100 mg/dL and HDL above 50 mg/dL are at no risk. In this paper, we developed a mathematical model of the formation of a plaque, which includes the following key variables: LDL and HDL, free radicals and oxidized LDL, MMP and TIMP, cytockines: MCP-1, IFN-γ, IL-12 and PDGF, and cells: macrophages, foam cells, T cells and smooth muscle cells. The model is given by a system of partial differential equations with in evolving plaque. Simulations of the model show how the combination of the concentrations of LDL and HDL in the blood determine whether a plaque will grow or disappear. More precisely, we create a map, showing the risk of plaque development for any pair of values (LDL,HDL).


Introduction
Atherosclerosis, hardening of the arteries, is the leading cause of death in the United States, and worldwide.The disease triggers heart attack or stroke, with total annual death of 900,000 in the United States [1] and 13 million worldwide [2].
Atherosclerosis is a disease in which a plaque builds up inside the arteries.A plaque contains low density lipoprotein (LDL), macrophages, smooth muscle cells (SMCs), platelets, and debris.The plaque constricts the lumen of the blood vessel thereby increasing the shear force of blood flow [3,4].As the plaque continues to grow, the increased shear force may cause rupture of the plaque, possibly resulting in the formation of thrombus (blood clot) [3,5], ischemic stroke, and heart attack [3][4][5].
The process of plaque development begins with a lesion in the endothelial layer, allowing LDL, to move from the blood into the intima and becoming oxidized LDL (ox-LDL) by free radicals (FRs).FRs are oxidative agents continuously released by biochemical reactions within the body, including the intima [6][7][8].Endothelial cells, sensing the presence of ox-LDL, secrete monocyte chemoattractant protein (MPC-1) [9,10], which triggers recruitment of monocytes into the intima [11].After entering the intima, monocytes differentiate into macrophages, which have an affinity for the ox-LDL [12][13][14].The ingestion of large amounts of ox-LDL transforms the fatty macrophages into foam cells [12,15].Foam cells secrete chemokines which attract more macrophages [10,12,13].SMCs from the media move into the intima by chemotactic forces due to MCP-1 [9,10], and platelet-derived growth factor (PDGF) [10,16], as well as by haptotaxis by the extracellular matrix (ECM).PDGF is secreted by macrophages, foam cells and SMCs [16,17].ECM is remodeled by matrix metalloproteinase (MMP) produced by a variety of cell types including SMCs [18], and is inhibited by tissue inhibitor of metalloproteinase (TIMP) produced by macrophages and SMCs [19].Interleukin IL-12, secreted by macrophages and foam cells [10,12,20], contribute to the growth of a plaque by activating T cells [9,20,21].Indeed, the activated T cells secrete interferon IFN-c, which in turn activates macrophage in the intima [13,21,22].At the same time that LDL enters the intima, high density lipoprotein (HDL) also enters into the intima, and becomes oxidized by free radicals [7,8].However, oxidized HDL (ox-HDL) is not ingested by macrophages.HDL helps prevent atherosclerosis by removing cholesterol from foam cells, and by the limiting inflammatory processes that underline atherosclerosis [23].Furthermore, HDL takes up free radicals that are otherwise available to LDL.Some of the key players in the atherosclerosis process are shown in Fig. 1.
It has long been recognized that the cholesterol concentrations in the blood are indicators of the probability that a plaque will develop: higher LDL and lower HDL concentrations indicate a higher probability of plaque development.Public health guidelines in the U.S. specify what levels of LDL are low risk and what levels are high risk; they also specify what levels of HDL are poor and what levels are near ideal [24,25].However, what is more relevant is to specify the risk associated with combined levels of LDL and HDL, and this is what the present paper addresses.A schematic of the network of atherosclerosis is given in Fig. 2. In this paper, we developed a mathematical model of plaque formation by a system of partial differential equations based on Fig. 2. The aim of the model is to determine the risk of plaque formation for combined levels of LDL and HDL.In particular, we created a ''risk-map'' for plaque development in the LDL-HDL coordinate plane, where the first quadrant of the plane was divided into regions of high risk, low risk and no risk.Anti-cholesterol drugs are aimed at lowering high levels of LDL, but some drugs are known to also increase the level of HDL [24].Hence such a risk-map may be important when evaluating the extend to which an anti-cholesterol drug can reduce the risk of atherosclerosis for particular individuals.

Mathematical model
In this paper, we present a mathematical model based on the network shown in Fig. 2. The model includes the variables listed in Table 1.We assume that all cells are moving with a common velocity u; the velocity is the result of movement of macrophages, T cells and SMCs into the intima.We also assume that all species are diffusing with appropriate diffusion coefficients.The equation for each species of cells X has a form where the expression on the left-hand side includes advection and diffusion, and F X accounts for various growth factors, bio-chemical reactions, chemotaxis and haptotaxis.The equation for the chemical species are the same but without the advection term.Fig. 3 shows a 2D cross section of a blood vessel with plaque V, and a planar cross section of a plaque in the direction along a blood vessel.
Equations for lipoproteins [LDL (L), HDL (H), ox-LDL (L ox )] and free radical (r).The distribution of LDL, HDL, ox-LDL and free radicals in the intima are described using reactiondiffusion equations [7], where k L and k H are reaction rates of oxidization, and l LoxM is the reduction rate of ox-LDL due to ingestion by macrophages.Eqs.
(1) and (2) model the evolution of LDL and HDL concentrations.
It is assumed that LDL and HDL are lost by reaction of oxidation with free radicals.Equation (3) models the production of ox-LDL due to LDL oxidation by reaction with the radicals (first term on right-hand side) and a reduction of ox-LDL through ingestion by macrophages (second term on right-hand side).Equation ( 4) models the evolution of free radicals concentration with baseline growth r 0 .Equation for macrophages (M).The evolution of macrophage density is modeled by Here the first term on right-hand side accounts for recruitment of macrophages by MCP-1 [9], and the second term accounts for the activation of macrophages by IFN-c [13,21,22].Equation for MCP-1 (P).The MCP-1 equation is given by where the first term on the right-hand side is the production of MCP-1 by endothelial cells, whose density is assumed to be constant, under the influence of ox-LDL [9].

Equation for T cells (T).
The density of T cells, which are primarily CD4 + T cells [20], satisfies the equation In this equation, we assume that T cells are activated by IL-12 in conjunction with MHC-II (major histocompatibility complex, class II).Actually, T cells are also activated by IL-1 and IL-6 produced by macrophages and SMCs [10,12,13].However, because of lack of experimental data, we do not include the IL-1 and IL-6 explicitly but instead consider their effect implicitly in estimating where the first term on right-hand side represents production of I c by T cells [21].
The first two terms on right-hand side account for chemotaxis by MCP-1 [9,10], and PDGF [10,16], and the last term accounts for haptotaxis by ECM.
The first term of right-hand side is the production of I 12 by macrophages enhanced by I c and resisted by HDL [23].The production of I 12 by macrophages is resisted by I 10 (which, for simplicity, is accounted by the factor 1 K M zM ) [10,12].The second term represents the production of I 12 by foam cells [20].
Equations for PDGF (G), MMP (Q) and TIMP (Q r ).We have the following sets of reaction diffusion equations for the chemokines (G, Q and Q r ): In equation ( 11), PDGF is produced by macrophages, foam cells, and SMCs [16,17].In equation ( 12), MMP is secreted by SMCs [18] (first term on right-hand side), and is lost by binding with TIMP (second term).In equation ( 13), TIMP is produced by SMCs and macrophages [19].
Equation for foam cells (F).Macrophages that have ingested a large amount of ox-LDL become foam cells [7,12,15], so we have Equations for ECM (r) and pressure (s).We assume that the intima has the constituency of a porous medium.Then, by Darcy's law, the velocity u of the cells is given by where s is the pressure.We also assume that the total density of all the cells plus the concentration of r is constant.This constant should be smaller than the average density of a plaque, 1.2260.03g/cm 3 [26], because plaques contain some debris, which are not included in our model.We take the constant to be 1 g/cm 3 , i.e.,

MzF zTzSzr~1: ð16Þ
We assume that all cells are approximately of the same volume and surface area, so that the diffusion coefficients of the all cells have the same coefficient, D. By adding Eqs. ( 5), ( 7), ( 9) and ( 14), we get where Eq. ( 17) gives a relation between r and s.We next derive an equation for r.The ECM is degraded by MMP [27], and is remodeled by macrophages and SMCs [18,27].For simplicity, we take the remodeling rate to be a constant, l r , as in [28].Then the equation of the density of ECM is given by Since u~{+s, this equation can be written in the form Adding this equation and ( 17), we get our final equation for s: The equation for r can then be written as Lr Lt {DrDrz(wzy)r{+s : +r~y: ð22Þ

Boundary conditions
For simplicity, we consider only 2-dimensional plaques as in Fig. 3. Then the boundary of the plaque consists of i) C M , in contact with media; ii) a free boundary C I , inside the lumen, and iii) two more vertical boundaries C L and C R of the intima in the case of Fig. 3(B).
Boundary conditions on C I .We assume flux boundary conditions of the form for X = L, H, M, T, and non-flux boundary conditions for all other variables, where Y ~Lox , r, P, I c , S, I, I 12 , G, Q, Q r , F. The boundary values for r are determined by Eq. ( 16).The coefficient a X is a constant except for M, and a M (L ox ,H)~ã a M Lox 1zH , since ox-LDL attracts monocytes [11], while HDL limits the inflammation process [23].Note that L 0 and H 0 are the LDL and HDL concentrations in the blood, so we shall be interested to see how these concentrations determine whether a small plaque will grow or shrink.
As in [29][30][31], we assume that the free boundary C I is held together by cell-to-cell adhesion forces so that where k is the mean curvature of the surface C I .(If C I is circular, then k is the reciprocal of the radius) Furthermore, the continuity condition u : n~V n , where n is the outward normal and V n is the velocity of the free boundary C I in the direction n, yields the relation Boundary conditions on C M .We assume non-flux boundary conditions for all variables except r and S on C M : For S, we have where S 0 is SMCs density in the media, and for simplicity we take a S ~aS (P,G)~ã a S PzG P 0 zG 0 since MCP-1 and PDGF attract SMCs from the media.As in the case of C I , the boundary values of r are determined by Eq. ( 16).Boundary conditions on C L and C R .We assume the periodic boundary conditions on C L and C R .

Parameter estimation
Table 2 lists the range of molecular weights of proteins and Table 3 lists their range of concentration.In the second columns in Tables 2 and 3, we indicate the (intermediate) values used in the simulations.The Tables 2 and 3 are used to estimate some of the model parameters.A summary of all the model parameters is given in Tables 4 and 5.
Reaction rates.To estimate some of the parameters in the equations for proteins, we shall use the concept of ''accessible surface area'' [32,33] of a protein p, or briefly ''area,'' A p , which is roughly the minimum surface area of the smooth shapes containing the protein.It was estimated in [34], that A LDL ~1:02|10 {11 cm 2 , and A HDL ~1:34|10 {12 cm 2 , so that their ratio is Accordingly, the corresponding reaction rates of the oxidation, k L and k H , are related by k L ~7:6k H .Moreover, the reaction rate of oxidation of LDL by free radicals is k L ~10:37 gcm {3 day 21 [7,34,35], so that k H ~1:36 g cm 23 day 21 .Diffusion coefficients.We assume that the diffusion coefficients of all the cells are the same, and take them to be 8:64|10 {7 cm 2 day 21 [28,36].In order to estimate the diffusion coefficients of the various proteins, we assume that the diffusion coefficient of protein p, D p , is proportional to its area A p , i.e., D p ~KA p , where we take K to be the same for all small molecules.For glucose, which is a monomeric globular protein, A p can be computed in terms of the molecular weight M p , by the formula A p ~11:1M 2=3 p [33], and M p ~180 dalton [28] The diffusion coefficient of MMP is 4:32|10 {2 cm 2 day 21 [38].
T cell 1610 23  Range of CD4 + T cells in healthy normal adult of 500,000 to 1,500,000 cells per ml [49].
doi:10.1371/journal.pone.0090497.t003PDGF is produced by SMCs, and likely also by endothelial cells and macrophages [41].In wound healing, macrophages produce PDGF at rate of 5.76 day 21 [42].Since the plaque formation is a much slower process, we take this rate l GM to be much smaller, i.e., l GM ~0:1 day 21 .Since SMCs production rate of PDGF is higher than that by macrophages [41], we take l GS ~0:5 day 21 .
The production rate of MMP by tumor cells was estimated in [28] to be 6|10 {3 day 21 .We assume that SMCs produce MMP at a much lower rate, namely, l QS ~3|10 {4 day 21 .Since SMCs produce MMP to enable them move into the intima by haptotaxis, we assume that they produce TIMP at a lower rate than MMP, and take l QrS ~1 10 l QS ~3|10 {5 day 21 .As macrophages produce most of the TIMP [43], we take the production rate of TIMP by macrophages to be l QrM ~2l QrS ~6|10 {5 day 21 .
We assume that the production rate of MCP-1 by endothelial cells, l PE , is twice that of d P P 0 , where P 0 is the concentration of MCP in the blood, which is equal to 1|10 {9 gcm {3 [44].We assume that d F ~1 2 d M ~0:03 day 21 [39], and that, in Eq. ( 14), l FM ½M&4d F ½F and ½F &½M, so that l FM ~4d F ~0:12 day 21 .
By [45], l I12M ~3|10 {7 day 21 .We assume that foam cells have lower production rates of I 12 and PDGF than macrophages, and take l I12F and l GF to be one third of the values of l I12M and l GM , respectively, so that l I12F ~1|10 {7 gcm {3 day 21 , and l GF ~0:033 day 21 .
Degradation rates.The degradation rate of MMP is d Q ~4:32 day 21 [37].Since TIMP has a short half life compared to MMP [46], we take its degradation rate to be d Qr ~5d Q ~21:6 day 21 .
In [47], the binding rate of MMP and TIMP is reported to be 3|10 5 M {1 s {1 , where 1M~the mass per mole, and the molecular weights of MMP and TIMP are 52 kda and 25 kda, respectively [48].Accordingly, we derive the binding rate per Molar per second (by same formula as in [44]), where N A is called the Avogadro number, and is the number of molecular per dm 3 .N A ~6:02|10 23 , and 1:66|10 {24 g is the mass of a proton for atomic mass unit.

Numerical methods
Finite element implementation.In order to illustrate our numerical method, we consider the following diffusion equation with Robin boundary conditions: where c+(uX ) is an advection term, and either c~1 or c~0 (no advection), and u~{+s.Multiplying the differential equation by an arbitrary function v, and performing integration by parts using the boundary conditions, we get This is an equivalent formulation of the system (27), which is better suited for simulation.
Similarly, Eq. ( 21) for s has the equivalent form: for an arbitrary function v. Galerkin discretization.The standard Galerkin discretization method uses finite dimensional subspaces V h to approximate the solution X.Let fv i g N i~1 be a basis of V h , where N is the number of nodes within the triangulation K. Let X n h [V h denote the numerical approximation of X at time t~ndt, where dt is the time step, X n h is written as a sum for coefficient c n j to be determined.
Recalling (30), we can rewrite the system (32) as a linear system of equations where c n is the vector of (c n j ), and the coefficient matrix A~(A ij ) and the right-hand side b~(b i ) are defined by Similarly, setting s n h ~PN j~1 d n j v j , where s n h is a numerical approximation of s at time ndt, Eq. ( 29) can be written as follows where d n is the vector (d n j ), B is a matrix (B ij ), B ij ÐV +v i : +v j dx, and e j ~D Ð V (+r : +v j )dxz Outline of the procedure.Suppose the domain V(t) has polygonal boundaries C I (t) and C O (t).Then we can cover the closure V(t) of V(t) by a regular triangulation T of triangles, i.e., V(t)~S T[T T where each T is a closed triangle.The triangular mesh, which is a basic thing that Finite Elements requires, is generated by distmesh [50], which is a mesh generation tool implemented in MATLAB, and our algorithm is outlined in Algorithm S1.For the detailed implementations, such as: construct basis functions over the triangulation, assemble the stiffness matrix, etc, see references [51,52].

Results
Numerical simulation is initialized by a small formed plaque.(see Figs. 4,5,6).Five combined levels of LDL and HDL (L 0 and H 0 ) are tested for 300 days:   Fig. 4 shows the growth of the plaque in case (a), Fig. 5 shows the shrinkage of the plaque in case (c), and Fig. 6 shows almost no plaque in case (e).In Fig. 7, the weight of the plaque, the summation of total cells, namely, Ð V (MzF zSzT)dV, is plotted for these five scenarios of combined levels of LDL and HDL.Similarly to Fig. 7, we show in supporting information files how the populations of macrophages, SMCs, foam cells and T cells, as well as the concentration of ox-LDL, IFN-c and IL-12, vary for different levels of LDL and HDL shown in Figs.S1, S2, S3, S4, S5, S6, S7.
Fig. 8 shows a risk-map of plaque development.To create the risk-map, we divided the LDL axis by 121 equidistant points, i.e, L i ~70zi (i~0, Á Á Á ,120), and divided the HDL axis by 21 equidistant points, i.e., H j ~40zj (j~0, Á Á Á ,20).For each pair (L i ,H j ), we computed the weight of the plaque, W (L i ,H j ) after 100 days on the domain shown in Fig. 3 (B), and formed the risk matrix where W 0 is the initial weight of the plaque.The vertical axis on the right of Fig. 8 shows the legend of the percentage of plaque growth or shrinkage.Accordingly, we divided the LDL-HDL plane into three regions: region I predicts high risk of plaque development, region III predicts no risk, and the intermediate region II predicts low risk.

Sensitivity analysis
In order to support the robustness of the simulation results, we ran sensitivity analysis on parameters which appear in the differential equations and in the boundary conditions.The parameters chosen are those whose baseline was somewhat crudely estimated while at the same time they seem to play an important role in the development of the plaque.Specifically, we chose all the 15 production rate parameters from the third box of Table 4, all the 5 influx rate parameters from the third box of  Table 5, and L 0 , H 0 .We list all these parameters with their range, baseline and unit in Table 6.
Following the sensitivity analysis method described in [53], we performed Latin hypercube sampling and generated 100 samples to calculate the partial rank correlation coefficients (PRCC) and pvalues with respect to the weight of the plaque after 300 days.The PRCCs are shown in Fig. 9, and all the p-values (not shown here) are less than 0.01.A positive PRCC (i.e., positive correlation) means that an increase in the parameter value will increase the weight of the plaque while a negative PRCC (i.e., negative correlation) means increase in the parameter will decrease the weight of the plaque.We note that l QS is positively correlated, as it should be.Indeed, if l QS is increased then MMP (Q) is increased (Eq.( 12)) so that ECM (r) is decreased (Eq.( 19)) and hence the plaque weight Ð (MzF zSzT)dV is increased (Eq.( 16)).As another example, note that l LoxM is negatively correlated.Indeed, if l LoxM is increased then L ox is decreased (Eq.( 3)), and a M in the boundary condition will decrease, leading to smaller M, and then to smaller T and F. Similar explanation can be given to the other parameters.
The most significant positively correlated parameters are L 0 and its influx rate a L .This is not surprising since LDL initiates the plaque formation.The most significant negatively correlated parameters are H 0 and its influx rate a H . Indeed, since HDL reduces the availability of free radicals, it plays an important negative role in plaque formation.

Discussion
Atherosclerosis is a disease in which a plaque builds up inside an artery.The process of plaque formation begins when, as a result of a lesion in the artery, cholesterols LDL and HDL enter the intima, and LDL becomes oxidized by free radicals.Upon sensing ox-LDL, endothelial cells secrete MCP-1 which attracts monocytes from the blood.As monocytes enter to the intima, they differentiate into macrophages that ingest the ox-LDL and become foam cells.Foam cells attract more macrophages, followed by T cells from the blood, and SMCs from the media.HDL reduces the available free radicals, as well as inflammation within the evolving plaque, thus HDL acts to block plaque growth.
Public health guidelines in the U.S. specify that LDL level of 100-129 mg/dL is near ideal, 130-159 mg/dL is borderline high, and 160-189 mg/dL is very high, whereas concentration of HDL above 60 mg/dL is best, and below 40 mg/dL for men or below 50 mg/dL for women is poor [25].An important question is how to evaluate the risk of atherosclerosis for a pair of LDL and HDL taken together.This question is addressed in the present paper.We built a mathematical model of plaque development by a system of partial differential equations.The model includes two parameters: L 0 , the level of LDL in the blood, and H 0 , the the level of HDL in the blood.
The model can simulate the evolution of a small plaque for any pair of values of (L 0 ,H 0 ).In Figs. 4, 5, 6, we simulated the plaque evolution over a period of 300 days.For example, one extreme case of L 0 ~190 mg/dL, H 0 ~40 mg/dL, the plaque doubled after 300 days; in another extreme case of L 0 ~70 mg/dL, H 0 ~60 mg/dL, the plaque disappeared after 300 days.We created a risk-map by taking sampling points of LDL and HDL values, and computing the weight of the plaque for each pair (L 0 ,H 0 ) after 100 days.The map shown in Fig. 8, indicates the percentage of plaque growth or shrinkage for any such pair.We accordingly divided the (LDL,HDL) quadrant into three regions: high risk, low risk, and non risk.
The need to consider the ratio of LDL/HDL in predicting coronary heart disease was suggested in a case study by [54].The American Heart Association considers the ratio (Tc=HDL)w5 to  indicate high risk of heart disease, and the ratio (Tc=HDL)v3:5 to be risk free [55], where Tc denotes the total cholesterol, which is calculated by the formula Tr (Tr~triglycerides): Table 7 shows the National Cholesterol Education Program (NCEP) guidelines associated with plaque buildup [56].Accordingly, for the five points (a)-(e) in Results we have: Some anti-cholesterol drugs, such as statins, lower LDL and at the same time also increase the HDL [24].It is important to know which drugs can best achieve the desired risk-free balance between LDL and HDL, that is, bring the individual's (L 0 ,H 0 ) into the no risk (or low risk) region.By focusing not on just reducing L 0 or on just increasing H 0 , but on moving the combined (L 0 ,H 0 ) to the no risk (or low risk) region in the shortest medically feasible path, we believe one could choose a more personalized medicine from those currently available, which will reduce the risk of atherosclerosis with the lowest amount of doze, thereby also possibly reducing potential negative side effects.
To illustrate this approach, we note that for some drugs the ratio of decrease in LDL to increase in HDL is already known.For example, this ratio is 1/3 for the new experimental drug Evacetrapib.Some anti-cholesterol drugs only decrease LDL (e.g.Colestid) while others only increase HDL (e.g.Lofibra).We can the represent effect of such drugs by unit vectors in the (L 0 ,H 0 ) plane: for example, Colestid r, Lofibra q, and Evacetrapid 9.In Fig. 10.we consider three individuals, A, B, and C in the high risk region.In order to move them to the low risk region with the minimum amount of medication (side effects are ignored), the individual should choose the drug for which the line segment from the individual initial position to the low risk region is the shortest (We assume that the amount of drug is proportional to the length of the line segment).Thus the best drug for A is the one that primarily increases HDL.Similarly, C will do better with a drug that primarily decrease LDL, and B should use a drug with appropriate ratio of decreasing LDL to increasing HDL.
Some work has been done on antioxidant therapy for reducing the risk of atherosclerosis, but so far it has had limited success in preventing cardiovascular diseases [57][58][59].A review of studies in which antioxidant gene therapy has been successfully used is given in [60].Our model could account for antioxidative medication once we gain a good understanding of how such medication affects the source of free radicals, r 0 , in Eq. (4).Some of the parameters in the differential equations in our model had to be rather crudely estimated, since no data were available, while others may slightly vary depending on the individual.As more data become available, parameter values may be further refined.Our model uses only the values of LDL and HDL as biomarkers.It will be interesting in the future to incorporate also triglycerides into the risk-map.Future work should also explore how other risk factors, such as high blood pressure, smoking and diabetes affect the risk-map.
We did not include in this paper the circulation ox-LDL in the blood, which is elevated only in patients with advanced atherosclerosis [61,62].Our model could be extended to include this additional biomarker, but at present there is not enough data on how the level of ox-LDL in the blood correlates to a specific advanced state of the disease.

Figure 1 .
Figure 1.Atherosclerosis schematics: the presence of ox-LDL in the intima causes monocytes to migrate from the lumen into the intima.Monocytes differentiate into macrophages which endocytose ox-LDL and become foam cells.SMCs are attracted from the media into intima by chemotaxis and haptotaxis.Cytokines released by macrophages, foam cells and SMCs activate T cells.T cells enhance activation of macrophages.HDL helps prevent atherosclerosis.doi:10.1371/journal.pone.0090497.g001

Figure 2 .
Figure 2. Schematic network of atherosclerosis.LDL and HDL are oxidized by free radicals, and become ox-LDL and ox-HDL respectively.Ox-LDL recruits macrophages to intima.By ingesting ox-LDL, macrophages are transformed to foam cells.SMCs are attracted into the intima by MCP-1 (secreted by endothelial cells) and PDGF (secreted by macrophages and foam cells).Macrophages, foam cells and SMCs secrete IL-12, which activates T cells.IFN-c secreted by T cells enhance the activity of macrophages which contributes the plaque built-up.doi:10.1371/journal.pone.0090497.g002

D
the parameter l TI 12 .For simplicity, we include the antiinflammatory effect of IL-10 produced by macrophages only implicitly, by the factor 1=(K M zM).Equation for IFN-c (I c ).The dynamics of IFN-c concentration is modeled by LI c Lt {D Ic DI c ~lIcT T |fflffl{zfflffl} production {d Ic I c |fflfflffl{zfflfflffl} degradation , LQ r Lt {D Qr DQ r ~lQrS Szl QrM M |fflfflfflfflfflfflfflfflfflfflfflffl{zfflfflfflfflfflfflfflfflfflfflfflffl} production {d QrQ QQ r |fflfflfflfflfflfflffl ffl{zfflfflfflfflfflfflffl ffl} depletion {d Qr Q r |fflfflfflffl ffl{zfflfflfflffl ffl} degradation

Figure 3 .
Figure 3. Two 2D cross sections of a plaque.C M is the boundary of the intima in contact with the media, and C I is the boundary of the intima in contact with the lumen.In (B) C L and C R are parts of the intima.doi:10.1371/journal.pone.0090497.g003

Figure 4 .
Figure 4. Simulations for the atherosclerosis model of 300 days after an initial plaque is formed with H 0 = 40 mg/dL and L 0 = 190 mg/dL.(A: Cross sections of a blood vessel, B:Cross sections along the blood vessel).doi:10.1371/journal.pone.0090497.g004

Figure 5 .
Figure 5. Simulations for the atherosclerosis model of 300 days after an initial plaque is formed with H 0 = 50 mg/dL and L 0 = 130 mg/dL.(A: Cross section of a blood vessel; B: Cross section along the blood vessel).doi:10.1371/journal.pone.0090497.g005

Figure 6 .Figure 7 .Figure 8 .
Figure 6.Simulations for the atherosclerosis model of 300 days after an initial plaque is formed with H 0 = 60 mg/dL and L 0 = 70 mg/dL.(A: Cross section of a blood vessel; B: Cross section along the blood vessel).doi:10.1371/journal.pone.0090497.g006 if Tr is above normal; (c) 3:5v Tc HDL v5 if Tr is not very high; (d) Tc HDL v3:5 if Tr is normal; (e) Tc HDL v2:5 if Tr is normal;.According to the NCEP guidelines, (a) and (b) should be in the high risk region; (c) in the low risk region; and (d), (e) in the no risk region, as indeed they are according placed in the risk map in Fig. 8.