Theoretical prediction and validation of cell recovery rates in preparing platelet-rich plasma through a centrifugation

In the present study, we propose a theoretical framework to predict the recovery rates of platelets and white blood cells in the process of centrifugal separation of whole blood contained in a tube for the preparation of platelet-rich plasma. Compared to previous efforts to optimize or standardize the protocols of centrifugation, we try to further the physical background (i.e., based on the multiphase flow phenomena) of analysis to develop a universal approach that can be applied to widely different conditions. That is, one-dimensional quasi-linear partial differential equation to describe the centrifugal sedimentation of dispersed phase (red and white blood cells) in continuous phase (plasma) is derived based on the kinematic-wave theory. With the information of whole blood volume and tube geometry considered, it is possible to determine the positions of interfaces between supernatant/suspension and suspension/sediment, i.e., the particle concentration gradient in a tube, for a wide range of centrifugation parameters (time and acceleration). While establishing a theory to predict the recovery rates of the platelet and white blood cell from the pre-determined interface positions, we also propose a new correlation model between the recovery rates of plasma and platelets, which is found to be a function of the whole blood volume, centrifugal time and acceleration, and tube geometry. The present predictions for optimal condition show good agreements with available human clinical data, obtained from different conditions, indicating the universal applicability of our method. Furthermore, the dependence of recovery rates on centrifugal conditions reveals that there exist a different critical acceleration and time for the maximum recovery rate of platelets and white blood cells, respectively. The other parameters such as hematocrit, whole blood volume and tube geometry are also found to strongly affect the maximum recovery rates of blood cells, and finally, as a strategy for increasing the efficiency, we suggest to dilute the whole blood, increase the whole blood volume with a tube geometry fixed.

Introduction For a few past decades, we have seen the increasing interest and advances in clinical applications of platelet-rich plasma (PRP) to various fields of plastic surgery, dentistry, orthopedics, sports medicine and so on [1][2][3][4][5]. These, in general, require the treatment of chronic wounds and/or muscle injuries, which can be greatly benefited from the positive potentials of PRP in the tissue healing and bone regeneration [4,[6][7][8][9][10]. In order to fully utilize the functionalities of PRP that make these applications promising, it is required to maximize the concentration of platelets (and/or white blood cells, although the beneficial effects from the inclusion of leukocytes are still under a debate), which plays a critical role in releasing growth factors, cytokines and proteinases [8,9,[11][12][13]. Thus, many previous studies have tried to optimize (or standardize) the protocols to prepare the platelet-rich plasma so far [14][15][16][17], to name some, while various commercial products to produce PRP have been introduced to the market and tested as well [18][19][20]. Since the quality and functionality of PRP are strongly dependent on the protocol used for its preparation, however, the wide variations in the reported conditions to prepare PRP, such as centrifugal acceleration and time, amount of volume of blood, and the type of anticoagulant platelet agonist, make it very difficult to compare the subsequent results fairly. Thus, a systematic relation between the preparation condition and the concentration of platelets (and/or white blood cells) is sorely required to clarify the clinical benefits (biological effects) of PRP. This issue is considered to be much more significant from the fact that the recovery rate of platelets from the commercial automated system, typically of high cost, is relatively lower (about 40-60%) than expected.
In general, a PRP preparation involves sequential steps of blood collection, centrifugation to separate and recover the platelets, and activation of the platelets. The centrifugation step, which is the main interest of the present study, consists of the first stage to separate red blood cells (RBCs) and the second one to concentrate platelets [15][16][17]21]. As shown in Fig 1, a whole blood (WB) is initially collected in a tube (with anticoagulants) and the first centrifugation is carried out at a constant speed to separate the RBCs from the whole blood. After this process, the WB is separated into three layers: an upper layer that is mostly occupied with plasma and platelets, an intermediate thin layer including a small amount of platelets and white blood cells (WBCs), and a bottom layer where most of the RBCs is packed. A typical PRP is usually obtained by transferring the upper and intermediate layers to an empty tube. The second centrifugation is then applied to concentrate the platelets and WBCs in the PRP. Thus, it is evident that the final properties of the PRP and efficiency of the whole procedure are strongly subject to the centrifugal conditions such as time (t c ) and speed (or centrifugal acceleration, a c ) [16,22], not to mention the influence on the maximized recovery of intact platelets (and/or WBCs).
As we have mentioned above, many previous studies have tried to find out the optimal condition of centrifugation for maximizing the recovery of platelets during PRP preparation [14][15][16][21][22][23][24]. For example, Slichter & Harker (1976) [22] showed that the single-step centrifugation of WB at a centrifugal acceleration of a c = 1,000g (g is the gravitational acceleration) for t c = 9 minutes provides the maximum recovery rate (up to 89%) of platelets. It was also reported that the shorter centrifugal time could be achieved by applying higher acceleration while the platelet viability is reduced significantly when the acceleration becomes greater than 3,000g, being applied for 20 minutes. Kahn et al. (1976) [14] explored the most efficient set of centrifugal conditions at relatively higher accelerations (a c = 1,614g-3,731g). They showed that the recovery rate of platelets does not appreciably increase when the spin time is longer than 8 minutes and the applied acceleration is larger than 2,324g. Considering the loss of platelet integrity together, thus several studies [15,16,21,23] recently reported optimal centrifugation conditions at relatively smaller accelerations (a c < 1,000g) to achieve the maximum recovery of platelets from the first stage of spin. Unfortunately, it is found that the results from these studies are inconsistent with each other and show a wide scattering in the optimal centrifugal conditions and the achievable recovery of platelets. Araki et al. (2012) [15], for example, obtained 80%-90% maximum recovery rate of platelets with a c = 230g-270g and t c = 10 minutes, while Amable et al. (2013) [23] reported about 87.7% platelet recovery rate with a c = 300g and t c = 5 minutes. Jo et al. (2013) [16] achieved a quite high platelet recovery rate of 92% at a c = 900g and t c = 5 minutes, and Perez et al. (2013) [21] obtained about 80% maximum recovery rate of platelets with a c = 100g and t c = 10 minutes. On the other hand, it is also necessary to solidify our understanding on the dependence of WBC recovery on centrifugal conditions, in addition to clearing the clinical role of WBCs [9]. Araki et al. (2012) [15] obtained about 25% recovery of WBCs and 70%-80% of platelets with a c = 70g and t c = 10 minutes. Perez et al. (2014) [24] showed the maximum recovery of WBCs (*10%) and platelets (*80%) with a c = 100g and t c = 10 minutes. It is noted that they all showed the recovery of WBCs is reduced to almost zero when the centrifugal acceleration is larger than about 800g. Now, it is clear that the widely scattered centrifugal conditions and resulting recovery of platelets and WBCs in PRP make it necessary to develop a universal protocol to optimize the PRP preparation method, which can be benefited significantly from a systematic analysis based on a physical theory underlying the separation process, i.e., multiphase flow phenomena. Considering the plasma and platelets (and/or WBCs) as continuous and dispersed phases, respectively, it is possible to utilize our knowledge in particle sedimentation. Perez et al. (2013) [21], sharing the same idea as ours, suggested a mathematical model to predict the recovery of platelets, based on the Stokes' law that describes a settling velocity of a single spherical particle in an unbounded environment. While the effect of backflow of the cell suspension has been included in their model, we think that more improvements are required to consider the Theoretical optimization of PRP preparation PLOS ONE | https://doi.org/10.1371/journal.pone.0187509 November 2, 2017 influences from particle-particle interactions and other operating conditions such as confined tube geometry and initial volume of whole blood.
Therefore, in the present study, we propose a new theoretical model to predict the recovery rates of the platelets and WBCs in PRP prepared by spinning whole blood contained in a tube. To solidify the physical background of the developed model, we apply the kinematic wave theory to estimate the volume fractions of each phase, i.e., plasma, platelets and WBCs, as a function of centrifugation time, acceleration, tube geometry, and the initial state of whole blood. To achieve this, we also perform a dimensional analysis to decide principal non-dimensional parameters that determine the recovery rate. The proposed model is applied to a wide range of parameters including volume of WB, hematocrit (volume fraction of RBC in WB), centrifugal time and acceleration, which is at the same time fully validated with available experimental (human clinical) data. Further discussions on the dependency of PRP composition on various parameters are given, and the results found in this study are expected to formulate practical guidelines for optimizing the PRP preparation.

Establishment of theoretical model
One-dimensional kinematic wave for the centrifugal sedimentation In the present study, we start from the idea that the packing of RBCs and WBCs at the bottom of the tube during the centrifugal separation of a whole blood can be considered as a centrifugal sedimentation of particles (dispersed phase) in a liquid (continuous phase). Especially, the one-dimensional kinematic wave theory is a very useful tool for our analysis, which is well known for its capability to determine gradient of particle concentration (i.e., interfaces between supernatant, suspension and sediment) in a tube [25][26][27][28][29][30]. Here, we consider the RBCs, WBCs and platelets as spherical solid particles of uniform diameter (d RBC , d WBC , and d PLT , respectively) and density (ρ RBC , ρ WBC , and ρ PLT , respectively), and the plasma as a liquid phase (density, ρ plas and dynamic viscosity, μ plas ). In Fig 2A, we have illustrated the present problem of spinning (at a constant angular velocity of ω) a tube that contains an incompressible mixture (i.e., initial whole blood, WB) of solid and liquid phases. Considering that a typical centrifugal acceleration (acting along the r-direction) applied for the blood separation is as large as O(10 2 -10 3 )g, the effect of gravitational field acting in z-direction is assumed to be negligible, and thus the particle movement along r-direction (i.e., one-dimensional motion) is considered to be dominant [26,28,29]. For the shape of tube bottom, we consider a flat geometry, simplifying that of a typical plain tube used in common [16]. As we will explain below, our prediction model can be applied to different tube geometry, i.e., the effect of tube geometry can be evaluated, and thus additional tube bottom shape is considered for the comparison ( Fig 2B). Now, we will constitute the mass and momentum conservation equations for the present problem. For a one-dimensional solid-particle suspension in a long tube, like the present case shown in Fig 2A, the unsteady continuity (mass conservation) equations for solid phase and mixture of liquid-solid phase are written as: and @j @r respectively, where α(t, r) is the volume fraction of the solid phase (i.e., particle concentration), j(t, r) and j s (t, r) denote the volume flux density [ms −1 ] of liquid-solid mixture and solid phase, respectively. Total volume flux density is then calculated as the sum of solid (j s ) and liquid (j l ) volume flux densities. Based on drift-flux model [31], the volume flux density of a solid phase in the reference frame moving at the volume flux density (j) of liquid-solid mixture, that is the drift flux j sl , is expressed as: where v s (t, r) denotes the velocity of a solid particle. Since the volume flux density j s is related to the concentration α such as v s = j s /α, then the Eq 3 can be rewritten as: Considering that the settling velocity of a particle should change along the r-direction, on the other hand, the momentum conservation equation for the solid particle motion relative to liquid phase can be reduced to a relationship between the drift flux (j sl ) and the settling function that is associated with the concentration α, as shown below [26,28]: Here, the function f bk (α) is defined as a Kynch batch flux density function [28] that describes the effect of existence of adjacent particles (i.e., effects from multiple particles are considered) on the settling (terminal) velocity of a single particle (u 1 ) through a stationary liquid, determined by the Stokes' law such as u 1 ¼ ðr s À r l Þd 2 s g=18m l (subscripts 's' and 'l' denote solid and liquid phase, respectively), where ρ and μ denote the density and viscosity, respectively, and d s is the size of the particle. That is, it accounts for the retarded settling due to the  Theoretical optimization of PRP preparation interactions between particles and a backflow from the bottom of the vessel. For the function f bk , we use the well-known Richardson-Zaki model [32], which has been widely used to similar problems to the present one [26,28,31]. The actual model is given as: where m is an index (or exponent) that is a function of the particle Reynolds number (Re p ) [31,32]. For a spherical solid particle, m = 5 is typically used when Re p is smaller than 0.2 [26]. For the present study, the particle Reynolds number (for RBCs and WBCs) is Re p * O(10 −2 -10 −1 ). This function also indicates that the settling is terminated at a limiting value of α = α max which is the maximum concentration of particles in the packed bed, i.e., in the bottom layer (see Fig  2A). In the present analysis, we use the maximum concentration value of α max = 0.8, considering the deforming property of the blood cells [27,29,33]. As introduced above, Perez et al. (2013) [21] also suggested a theoretical prediction of cell recovery, based on the modified settling velocity of a particle under centrifugation. They have corrected the settling velocity (u 1 ) of a single particle (i.e., RBC) to consider the effect of backflow as: where K (= 4.87) is a constant and α BL indicates the concentration of RBCs in the bottom layer, whose value was determined based on the conservation of RBC concentration between the WB and bottom layer in their study. Therefore, this approach has a limitation such that the value of constant K needs to be found for each condition. Furthermore, unlike the present model, the effect of particle interactions was not been explicitly included in this model. By substituting the Eqs 4 and 5 into the continuity equations (Eqs 1 and 2) while imposing j = 0 at r = R o (due to the closed bottom wall of the tube, see Fig 2A), we obtain This is further non-dimensionalized by introducing the dimensionless variables of r Ã = r/R o , f Ã bk ¼ f bk =u 1 and t Ã = tω 2 u 1 /g. Now, the final dimensionless form of our governing equation for the particle concentration is expressed as: which is a first-order quasi-linear partial differential equation that describes the centrifugal particle sedimentation. Instead of using the method of characteristics [26,30], we solve the governing Eq 9 numerically to determine the positions of the interfaces between supernatant/ suspension (I n − s ) and suspension/sediment (I s − d ) (Fig 2A). Actually, the analytical solution for Eq 9 is not in a completely explicit form and requires numerical integration as well, which is furthermore known to be highly dependent on the range of initial condition (particle concentration) [34]. Thus, being free from any constraints, we decided to get the solutions numerically for a wide parameter (initial condition) range. As a result, after spinning for certain time duration, it is possible to distinguish the concentration zones such as a clear liquid (α = 0), the retarded settling zone (0 < α α max ) and the packed bed (bottom layer) (α = α max ). To numerically solve the Eq 9, on the other hand, we use the modified upwind finite difference method [35] and the Engquist-Osher scheme [36] for extrapolating the flux density function. The details of the applied numerical algorithms are explained in elsewhere [35]. , obtained by solving Eq 9, which is represented by the iso-concentration lines. Also, the WBC concentration distribution (see below for the method to predict it) under the same centrifugal condition is shown together. Here, it was assumed that the liquidsolid mixture is composed of plasma and RBCs (or WBCs) only (the maximum and initial concentrations of RBCs are chosen to be α max = 0.65 and α o = 0.4, respectively), and the applied centrifugal acceleration is fixed as a c = 1,000g. The considered tube geometry, initially filled with the mixture of plasma and RBCs up to L, is also shown together in Fig 3. As we have explained above, it is now possible to determine the interfaces between regimes with different RBC concentrations. For example, after the tube is spun for 100 seconds, three distinct interfaces are detected depending on the variation in RBC concentration (Fig 3). In particular, for the case of monodisperse biosuspension like the present plasma-RBCs mixture in a tube, an upper layer (so-called PRP) consisting mostly of plasma can be determined in terms of the interface (I n − s ) between the region without any RBCs (pure plasma) and that with RBCs. Thus, the recovery rate of the plasma (E plas ) is calculated by [14,21]: where V UL (= (L − I n − s )A) is the volume of the upper layer (A: cross-sectional area of tube (= πD 2 /4) and L: height of mixture), and V WB and H e denote the volume of whole blood and hematocrit, respectively. While the Eq 10 was derived from the reasoning that RBCs take up about 99% of the total cells in WB [21], it is also clear that the upper layer, typically called as a PRP in practical centrifugal separations, contains platelets as well as plasma. This is because the platelets, whose nominal size and density is about 1.0 μm and 1,050 kg/m 3 , which are much smaller than those (8.0 μm and 1,125 kg/m 3 ) of RBC, would not migrate further to the tube bottom. Therefore, it would be possible to draw the relation between the recovery rate of platelets (E PLT ) and that of plasma (E plas ), which will be discussed in the next section. In the above formulation, we have approximated the blood cell as a rigid sphere to calculate the settling velocity (u 1 ) of a particle (Stokes' law). As is well known, this assumption has been used widely by previous studies on the centrifugal cell separation (or elutriation) [21,27,29,[37][38][39][40].
In these studies, they modified the viscosity, density or the relative velocity of particles to compensate the uncertainties that may arise due to assuming blood cell as a sphere. For the same purpose, actually, we have used the modified settling velocity model by using the flux density function [32] (Eq 6). Due to the difficulty in the direct evaluation the assumption, instead we compare the volume of serum (i.e., the upper layer in Fig 3), V UL = (L − I n − s )A, that is used to determine the recovery rate of plasma (Eq 10), obtained from experiment [16] and our prediction for the hematocrit range of H e = 0.37-0.4. As shown in Fig 4, the V UL 's both from experiment and prediction show a good agreement each other, indicating that the assumption spherical particle is reasonable. On the other hand, it was reported that the density of red blood cell (ρ RBC ) tends to vary during in vivo aging [41]. However, we have confirmed that the For the estimation of the recovery of white blood cells (WBCs), on the other hand, we start from the assumption that the particle-particle (i.e., RBCs-WBCs) interaction would be weak during the centrifugation, since the amount of WBCs is much smaller than that of RBCs. Then, the interface between pure plasma and WBC-dominant layer (I pW ), and that between WBC-dominant and RBC-dominant layers (I n − s ) is also determined by solving the governing equations (Eq 9) for WBCs and RBCs, separately. The exemplary locations of both interfaces that can be detected from the present approach are shown in Fig 3. It is again noted that the actual functionality of having WBCs in PRP is not the issue of present study. As a result, the recovery rate of WBC (E WBC ) is simply obtained by where c w is an empirical coefficient with an order magnitude of O(10 −3 -10 −2 ), introduced from previous studies of Sartory (1977) [42] and Berres et al. (2003) [43], and β o (= 0.01) is the typical initial volume fraction of WBC in WB [44].

Correlation model for platelet recovery rate
As we have explained in previous section, it is necessary and possible to draw a correlation between the recovery rates of platelets and plasma. Indeed, Brown (1989) [38] has claimed that the recovery rate of platelets is linearly proportional to that of plasma during centrifugal cell separation, i.e., E PLT /E plas = constant (' 1.0). Based on the theory of scales of measurement [45,46], on the other hand, the quantity E plas and E PLT are classified as the ratio scale-types, and the admissible relation between the independent quantity E plas and dependent quantity E PLT can be expressed as a product form of For conveniences, logarithmic form of Eq 12, a log(E plas ) + log(E PLT ) = log(b) is usually used, and a is a constant that is independent on any other quantities while b is dependent on other quantities. So the correlation model by Brown (1989) [38] corresponds to the case of a = −1 and b = 1 such that log(E PLT /E plas ) = 0, indicating that the ratio of recovery rate of platelets to plasma is independent on any other variables. While gathering and analyzing the actual clinical data [14-16, 21, 23] available in the literature, however, we found that the above relation is not always valid. Fig 5 shows the scattering of actual values of log(E PLT /E plas ) gathered from the above literatures together with that of Brown (1989) [38]. Considering the fact that the ranges of centrifugal time and acceleration, tube geometry, and the initial condition of tested blood are widely different between the referred studies, we think it is quite necessary to draw a more reliable (i.e., applicable to a more wide range of conditions) model to correlate the platelet recovery rate to that of plasma. To achieve this, we first perform a dimensional analysis, i.e., Buckingham Pi-theorem [47,48]. Based on our physical intuition and results from previous studies, it is reasonably deduced that the recovery rate of platelets (E PLT ) is determined by various variables as shown: Here, u 1 ¼ ðr s À r l Þd 2 s g=18m l and U o ¼ ðr s À r l Þd 2 s ðo 2 R o Þ=18m l are the terminal (settling) velocity of a single particle under the gravitational and centrifugal accelerations, respectively, and V t is the volume of the tube. The above functional form also indicates that log(E PLT /E plas ) should be expressed as Eq 14, rather than being a constant of 0 as suggested by Brown (1989) [38].
Here, the centrifugal time (t c ), cross-sectional area of the tube (A) and the initial volume of whole blood (AL) are considered as fundamental kinds of units [47] based on their physical importance in the relation. To constitute dimensionless groups (i.e., P's) with a physical dominance using Buckingham's Pi-theorem, first of all, we need to classify the variables in Eq 14 into primary and repeating ones, respectively. Since the major goal of present model is to optimize the centrifugal parameters and initial conditions in PRP preparation, it is reasonable to consider V t , U o and t c as variables of interest, while taking A, L and u 1 as repeating variables. Thus, from Pi-theorem, we can derive four P's as: where P 1 denotes the recovery rate ratio of platelets to plasma, P 2 is the volume fraction of WB relative to the tube capacity, P 3 shows the characteristic time scale, and P 4 explains the ratio of centrifugal acceleration to gravitational one (g), respectively. Consequently, a generalized correlation between Pi-groups based on a power function is expressed as: where the coefficients (c 1 and c 2 ) and exponents (e 2 , e 3 and e 4 ) are determined through the regression method with existing experimental data. In the present study, we used the experimental data of Araki et al.  [21], in which all the information necessary to complete the above relation was available (the validation will be discussed later). That is, as shown in Fig 6, we plotted the variations of P 1 with P 2 P 3 P 4 for the available experimental data from different conditions and they collpse into a single linear function of P 1 = −0.0122P 2 P 3 P 4 + 0.5128. As a result, the equation between Pigroups is established as: As shown, all the exponents are 1.0, which indicates that the logarithmic ratio of E PLT to E plas has a linear relationship with the dimensionless variables (P 2 , P 3 and P 4 ) when any two of them is fixed. Finally, the correlation model between E PLT and E plas is obtained as: Theoretical optimization of PRP preparation As noted, this correlation describes E PLT /E plas as a function of the whole blood volume, centrifugal time and acceleration, with which the platelet recovery rate (E PLT ) can be determined once the recovery of plasma (E plas ) is obtained by solving the Eq 9.
To further discuss the reliability of our model, we have compared the values of E PLT /E plas measured under different conditions [15,16,21] in Fig 7 with our predictions where each clinical environment has been fully incorporated. It is found that for each data set, the recovery rate ratio shows an exponential decay with increasing centrifugal acceleration, which agrees well with our model. Furthermore, the decaying rate strongly depends on the centrifugal time (t c ) and whole blood volume (V WB ). As the centrifugation time increases while fixing V WB , for example, the decaying rate becomes much faster; however, the increase of V WB (with a fixed t c ) would reduce the rate. This shows that we may obtain a high E PLT /E plas with a short centrifugation time or large blood volume with the same centrifugal acceleration, and there should be an optimal combination of controlled parameters to maximize the recovery rate of platelets. Later, we will discuss this condition in detail. Finally, we would like to emphasize that it was not possible to capture all these features with the previous correlation model by Brown (1989) [38], but our new model can successfully explain the physical details involved in the centrifugal separation of blood cells.

Determination of the parameter ranges
To cover various conditions in practice, a range of parameters such as volume of whole blood (WB), tube geometry (e.g., cross-sectional area), hematocrit (volume fraction of RBCs in WB) and centrifugal condition (centrifugal time and acceleration) are considered in the present  [16], respectively, while the typical tube geometry such as A = 120 mm 2 and V t = 15 mL is considered. The distance between the center of rotation and tube bottom is R o = 150 mm. To investigate the effect of tube geometry, we additionally consider the conical bottom shape as well (Fig 2B). We vary the hematocrit in the range of H e = 0.37-0.52, reflecting the difference between individuals (e.g., difference between male and female) and also the maximum concentration of RBCs due to packing at the tube bottom is limited by defining α max = 0.80, caused by the deformation of the blood cells [27,29,33]. This α max is defined such that the flux through the interface bordering the packed bed saturates to zero and for biosuspensions (such as RBCs) it has been set to be 0.8 [27,33,49]. Finally the ranges of centrifugal time and acceleration are determined based on the typical values used in clinical practices; so that the centrifugal time is varied as t c = 2-15 minutes and the centrifugal acceleration covers the range of a c = 100g-1,500g (g = 9.81 m/s 2 ) [14-16, 21-24, 50].

Comparison between prediction and experimental data: Validation
In this section, we will compare the predicted recovery rates of platelets and WBCs with the experimental data available in the literature. Fig 8 shows the variation of the recovery rates of platelets (E PLT ) with centrifugal acceleration (a c = ω 2 R o ) for three cases of V WB = 3.5, 7.0 and 9.0 mL, while the centrifugal time is fixed as t c = 10 minutes. Here, the compared experimental data were adopted from three different studies [15,16,21], indicating that they were obtained under different environments that were reflected in the prediction with the present model. Although almost all information to run our model were available, the exact value of hematocrit (H e ) was not specified and thus we used its typical range of H e = 0.37-0.52, of which the boundary is shown as dashed and solid lines in each plot in Fig 8. It is clearly shown that the present model predicts the dependency of E PLT on the centrifugal acceleration pretty well; in particular, the critical centrifugal acceleration (i.e., the spinning speed of centrifugation) for the maximum E PLT matches well with the experimental data. As a smaller volume of whole blood is used, the maximum E PLT achievable from the centrifugation is reduced, while the critical centrifugal acceleration decreases as well. In the following sections, we will discuss the effects of these principal variables in detail. On the other hand, the predicted E PLT shows a slight deviation from the experimental data for the centrifugal accelerations smaller than about 100g (see Fig 8B, for example). This is because the contribution of the gravitational acceleration would not be negligible as the centrifugal acceleration becomes smaller [26,30]. In formulating our model, we ignored the effect of gravitational force that acts along the direction perpendicular to the sedimentation direction (see Fig 2A).
In Fig 9, we also validated the accuracy of predicted recovery rates of WBCs (Eq 11). Unlike the case of platelet recovery, less experimental data were available for the inclusion of WBC's in PRP and single data set from Araki et al. (2012) [15] was compared in the present study. This is because the separation and recovery of WBCs during PRP preparation has not been paid attention much due to the controversy about the clinical role of WBCs in PRP [9]. In the figure, the solid lines denote the boundaries of the predicted E WBC for the range of coefficient c w (10 −3 -10 −2 ) [42,43]. It is found that our prediction captures the same trend of varying E WBC according to the centrifugal acceleration; in particular, the critical acceleration to achieve the maximum recovery of WBCs agrees well with the experimental data. Similar to the platelet recovery, the deviation between our prediction and experiment becomes larger at smaller centrifugal acceleration (< 100g). To our best knowledge, this is the first attempt to predict the recovery rate of WBCs in PRP preparation and the results are quite promising.

Effects of centrifugal conditions on the recovery of blood cells in PRP
Now, we discuss the effects of centrifugal conditions (i.e., time t c and acceleration a c ) on the recovery rates of platelets (E PLT ) and WBCs (E WBC ) based on the theoretical predictions developed in the present study. Shown in Fig 10 are the variations of the averaged (for the ranges of H e = 0.37-0.52 and c w = 10 −3 -10 −2 , respectively) recovery rates of platelets and WBCs with varying t c and a c . The maximum concentration of blood cells is set as α max = 0.8 and V WB = 9 mL is considered. When the spinning speed increases with a fixed time, both E PLT and E WBC increase sharply to approach the maximum value and then decreases relatively slowly. This dependency on centrifugal acceleration can be also found in the recent experimental studies https://doi.org/10.1371/journal.pone.0187509.g008 [15,16,24], but they considered a single value of t c while our predictions make it possible to perform a parametric study. As shown in Fig 10A, when t c is as small as 2 minutes, E PLT increases slowly and does not reach the maximum even with the acceleration as large as a c = 1,500g. This is because the centrifugal time is not long enough to separate the plasma and RBCs successfully. As t c increases, the slope of increasing E PLT becomes steeper and a larger recovery of platelets (over 80%) is possible at a critical a c that tends to get smaller as the centrifugation is performed longer. This information would provide us a good strategy to achieve the maximum platelet recovery by optimizing the condition for practical application. For example, the critical a c becomes as small as about 350g at t c = 12 minutes, which will be very helpful in maintaining the platelet integrity [14,22,24] as well as maximizing the platelet recovery. Interestingly, the maximum recovery rate of platelet achievable at each critical a c does not change much with varying t c , which may indicate that the maximum E PLT depends on the other conditions such as hematocrit, whole blood volume and tube geometry rather than the centrifugal time and acceleration (details will be discussed below). Fig 10B shows the effects of centrifugal conditions (t c and a c ) on the recovery rate of WBCs (E WBC ). The general trends of E WBC with t c and a c are similar to those of E PLT ; that is, with increasing a c , E WBC increases to the maximum at a critical a c , and then decreases back. The slope of E WBC change with a c and the critical a c show similar variations with those of E PLT , and the maximum E WBC is maintained constant despite varying centrifugal conditions. On the other hand, it is noted that the critical a c to achieve the maximum E WBC is smaller than that for the maximum E PLT , which agrees with the experimental data of Araki et al. (2012) [15]. This is because the time (velocity) scale of the platelets is much smaller than that of WBCs due to the large difference in their sizes and densities. Again, we think this is one of advantages of our theoretical model to capture this difference accurately. Theoretical optimization of PRP preparation After the feasibility of our model is confirmed, we will discuss more on the details of our model. First of all, to understand this discrepancy between the present model and Brown's assumption [38], it would be meaningful to collect the estimated recovery rates of platelets and plasma as a polar plot, as shown in Fig 11. The data set shown in the figure includes the cases of t c = 2, 5, 10, and 15 minutes with a c = 100g-1,500g. It is interesting to see that all the collected data collapse into a single parabolic curve. As we have introduced above, Brown (1989) [38] has proposed that the polar plot between E PLT and E plas has a linear correlation, which actually corresponds to the regime before the maximum E PLT is achieved (Fig 11). That is, in our analysis, E PLT is found to increase almost linearly with increasing E plas until it reaches the maximum. Thus, we may imagine that the Brown's result has been established based on the incomplete data set. After the maximum E PLT is achieved, it begins to decrease with further increase of a c ; on the other hand, E plas steadily increases. In spite of different centrifugation times, the collapse of the data into a single curve indicates that the maximum E PLT and the corresponding E plas are maintained the same. This suggests that there exists a critical E plas , determined by the volume of the upper layer, so-called PRP, in which the amount of recovered platelets is maximized. It is again shown that E plas increases with increasing t c , which would Theoretical optimization of PRP preparation reach the critical E plas at lower a c , indicating that the critical a c for the maximum E PLT decreases with increasing t c .

Effects of whole blood volume, hematocrit and tube geometry
The recovery rate of blood cells would be also strongly influenced by the conditions such as the volume of whole blood, hematocrit, and tube geometry. Fig 12 shows the variations of E PLT and E WBC with a c (t c is fixed as 10 minutes) while varying the volume of whole blood as V WB = 3.5-9.0 mL and hematocrit as H e = 0.37-0.52. First of all, the effect of V WB on the averaged (for the range of H e = 0.37-0.52) E PLT and E WBC is plotted in Fig 12A and 12B. When V WB is as small as 3.5 mL, E PLT reaches the maximum (*63%) at a very low centrifugal acceleration of 120g that is close to the optimal a c of 100g measured by Perez et al. (2013) [21] who used the same V WB . On the whole, with increasing V WB , the maximum achievable E PLT and corresponding critical a c tend to increase (Fig 12A). Thus, even considering the limited capacity of a tube (typically 15 mL in real applications), it is recommended to collect more blood to maximize the recovery rate of platelets. The recovery of WBCs shows a similar trend such that the maximum E WBC (and critical a c ) increases as V WB increases ( Fig 12B). As we have shown in Fig 10, the magnitude of E PLT is usually higher than E WBC , which has been measured in previous experimental studies [15,21], as well. This is because not only do WBCs have a very small Theoretical optimization of PRP preparation (< 0.01) initial volume fraction, but they also exist in a thin intermediate layer in PRP such that it is easy for them to be lost into the packing bed of RBCs. Finally, it should be again noted that the critical a c for the maximum E WBC is not the same as the one for the maximum E PLT , which is useful to choose the most efficient set of centrifugal conditions at a given clinical condition.
Shown in Fig 12C are the variations of E PLT with a c for three difference cases of hematocrit as H e = 0.37, 0.45 and 0.52. Again, these values of hematocrit were chosen to reflect the variations in real human data [52]. Other parameters are fixed as V WB = 9 mL and t c = 10 minutes.
As H e increases, the maximum E PLT decreases significantly while the corresponding critical a c increases. In this sense, it would be helpful to have a smaller volume fraction of RBCs in the whole blood to have more platelets in PRP, which explains the reason why previous studies have diluted the blood before centrifugation [17,21]. Similarly, the maximum E WBC and the corresponding critical a c are predicted to increase and decrease, respectively, as the initial volume fraction of RBCs decreases ( Fig 12D). As shown, the change rate of E WBC with a c near the maximum peak becomes much steeper as H e decreases, indicating that the E WBC is more sensitive to a c at lower H e . Therefore, we may imagine that controlling the separation and recovery of WBCs in PRP would be very difficult in the actual clinical environment.
Finally, to compare the effect of bottom shape of tube (flat bottom as a reference), we have additionally tested the conical bottom shape, which has been adopted from the shape of BD Falcon conical tube that is widely used in the actual clinical test. Typically, the height of conical part is about one fifth of the total length and the diameter of upper base has one quarter of the diameter of the tube. In the present study, the conical shape with the height of 0.2R o (R o , distance between the center of rotation and tube bottom) and the diameter of upper base of 0.25D (D, width of the tube) has been considered (see Fig 2B). Fig 13 shows the effect of tube geometry (i.e., flat and conical shapes at the bottom) on E PLT and E WBC , while we vary H e as 0.37 and 0.52 with V WB = 9.0 mL and t c = 10 minutes. In general, the overall dependency of Theoretical optimization of PRP preparation E PLT and E WBC on a c is not affected much by varying the cross-sectional area of tube bottom; however, the maximum achievable E PLT and E WBC (and corresponding critical a c ) have changed. Compared to the flat type, the conical type produces smaller maximum recovery of blood cells and the corresponding critical centrifugal acceleration becomes larger. For example, for the flat bottom geometry, E PLT reaches the maximum (*100%) at a c = 300g with H e = 0.37, but that of the conical tube has the maximum (*98%) at a c = 450g (Fig 13A). As a c increases toward the critical value, the predicted E PLT (and E WBC ) in the flat type is higher than that in the conical type. Interestingly, after the maximum recovery is achieved, this phenomenon is reversed at a c between 800g and 900g. To understand the reasons for this (especially focusing on E PLT ), we perform a scaling analysis based on the Eqs 10 and 18 that were used to determine the recovery rate of platelets. As shown below, it is possible to draw simple scaling relations for the change rates of E plas and E PLT /E plas on a c .
Once all variables for centrifugal separation of blood cells are pre-determined, the first relation shows that @(E plas )/@a c is a positive constant indicating the gathering rate of platelets in the PRP. On the other hand, the second relation, a decreasing function of a c in negative values, Theoretical optimization of PRP preparation denotes the rate of platelet loss to the packed bed regime. Thus, in Fig 14, we have plotted the variations of above two quantities (in absolute values) for the cases of H e = 0.52 shown in Fig  13A. Here, we can find a few interesting things. First of all, for both bottom shapes, it is shown that the two ratios cross each other at a c = 800g-900g (= a cr ) which is close to the centrifugal acceleration where E PLT 's from both cases are reversed (see Fig 13A). So, it is understood that the loss of E PLT would be faster at a c < a cr while more platelets will be gathered at a c > a cr . Second, the change rates of both E plas and E PLT /E plas with a c are higher for the flat type, compared to the conical one. Therefore, as shown in Fig 13A, the platelet recovery rate on flat type will be higher before the above two rates are balanced (i.e., a c < a cr ) while the conical one will prevail at a c > a cr .

Optimization of centrifugal conditions for PRP preparation
In previous sections, the effects of centrifugal conditions, whole blood volume, hematocrit, and tube geometry on recovery rates of platelets and WBCs have been discussed. Among these parameters, it is most frequently required to optimize the centrifugal conditions in real clinical applications. By mapping the predicted E PLT and E WBC together, it is possible to provide a simple guideline to optimize them. For example, Fig 15 shows the contours of averaged (for the range of H e = 0.37-0.52) E PLT and E WBC with a c and t c (whole blood volume is fixed as 9 mL and the flat tube bottom geometry was considered). It is clearly observed that there are two regions for large E WBC (! * 15%) and E PLT (! * 80%), and depending on the target it is feasible to select the appropriate centrifugal conditions. If one wishes to achieve the time efficiency, about 82% E PLT and 10% E WBC would be obtained at the centrifugal acceleration of 1,500g for a relatively short spinning time of t c = 2.5 minutes, or about 70% E PLT and 15% E WBC are obtained at a c = 1,200g and t c = 2.0 minutes. On the other hand, to achieve the enhanced platelet integrity, it is typically recommended to apply smaller centrifugal acceleration, which will require longer centrifugation time to maximize the recovery rates. For example, about 81% E PLT and 11% E WBC are achieved at a c = 240g for t c = 15 minutes, or about 63% E PLT and 16% E WBC are at a c = 140g for t c = 15 minutes. In general, it is understood that the trade-off between the time efficiency and platelet integrity should be considered in practical applications. In this case, the optimal set of centrifugal conditions consisting of the moderate centrifugal acceleration between 700g-800g and relatively short spinning time (4-5 minutes) would be selected to obtain 80% E PLT and 12% E WBC . Once the important environments are fully considered in the prediction, such as tube geometry, hematocrit and volume of whole blood, that may cause the large discrepancy between the optimal centrifugal conditions reported from previous studies, it is possible to generate this kind of map to select the spinning time and speed according to the specific cases. We would like to finally emphasize that this optimization was not possible in previous attempts.

Conclusion
In the present study, challenged by the large discrepancy between the recently proposed optimal centrifugal conditions for platelet-rich plasma (PRP) preparation, we developed a theoretical model, based on the knowledge in multiphase flow, to predict the recovery rates of the platelets and white blood cells in the process of centrifugal separation of the whole blood in a tube for the preparation of PRP. For a range of practical parameters, such as whole blood Theoretical optimization of PRP preparation volume, hematocrit and centrifugal conditions (i.e., time and acceleration), the predicted recovery rates show good agreements with available experimental data. The dependence of the recovery rate (platelets and WBCs) on centrifugal conditions is found such that there exist the critical time and acceleration to achieve the maximum recovery of platelets or WBCs. In addition, the effects of whole blood volume, hematocrit and tube geometry on the recovery rate have been investigated, confirming that the magnitudes of maximum recovery rates of platelets and WBCs are strongly affected depending on them. Furthermore, we have shown that the dilution of whole blood, increase of whole blood volume and extension of cross-sectional area of the tube are helpful to significantly increase the recovery rate. In view of the real clinical applications, our predictions can be used as a simple guideline to customize the centrifugal conditions according to the trade-off between practical issues such as the time efficiency and platelet integrity. We believe that the present results will explain the reason why the previously reported optimal conditions show wide scatters, and our model may satisfy the necessity of universal optimal condition to maximize the recovery of blood cells in PRP preparation, and will be provide a theoretical background to develop a specially customized (complex, for example) tube geometry with a enhanced capability of blood cell separation. Since our major goal was to predict the concetrnation of platlets (most important blood cell to determine the quality of PRP), we have used a relatively simple model to predict the WBC concentration. If the clinical role of WBCs in PRP is cleared, however, it would be also important to develop more rigorous model to estimate the WBC concentration in PRP preparation, as a future work.