The Effect of a Non-Homogenous Erythrocytes Properties Distribution in In-Silico Blood Flows

(1)

Bachelor Informatica

The Effect of a Non-Homogenous

Erythrocytes Properties

Dis-tribution in In-Silico Blood Flows

Sietse Molenaar

July 6, 2018

Inf

orma

tica

—

Universiteit

v

an

Ams

terd

am

(2)

Abstract

Blood is a vital part of the human body that performs many important functions. It is a non-Newtonian suspension yielded by deformable cells in a Newtonian plasma. The huge amount of tiny cells inside a blood suspension defines its behaviour and this makes modelling blood a complex task. For some biological phenomena looking at just the macroscopic scale of blood is not enough and observing what happens on a cellular level is of great importance. In this study an already existing blood simulation model is investigated further. This model can simulate thousands of cells in human blood on a microscopic timescale. The properties of these cells however are exact copies of each other and this research project examines what different distributions in these properties such as cell volume can be implemented and the effect this has on the blood flow properties. Distribution in cell size enables the ability to study the effects of conditions like iron deficiency and macrocytic anemia on both blood flow and transport properties.

(3)

Introduction

1.1 Context

A computing framework called HemoCell has been developed that simulates cells on a cellular level [1]. This allows researchers to investigate blood properties in different environments on a microscopic scale and with small time steps. The framework models not only RBC (red blood cells) but also platelets and white blood cells. The framework is written in C++.

To make the framework viable considering both development time and performance some sim-plifications have been made. For example RBC have properties like rigidness and size which are static values in the current HemoCell implementation. This means that when a simulation is done with thousands of RBC, all those RBC will have the same values for those properties. However, this does not match with physiological blood in which the properties are different for each cell [2].

The goal of this research project is to see how property distributions can be incorporated in the current version of HemoCell and what the effects of these changes are on the blood flow properties.

(6)

(7)

CHAPTER 2

Theoretical background

This chapter will give an overview of the relevant aspects of fluid dynamics and biomechan-ics.

2.1 Fluid dynamics

Just like other substances fluids are described by several properties. One of the most basic properties is viscosity. Viscosity is the measure of its resistance to flow under an applied force [3]. For example water is a fluid with a very low viscosity while oil has a slightly higher viscosity and blood has a viscosity somewhere in between [4]. Peanut butter and honey are examples of fluids with a very high viscosity.

Viscosity is defined using the ratio of the shear stress to the shear rate [4]. Shear forces are forces that push one part of a body in one direction and another part in another direction which causes deformation as seen in figure 2.1. The deformation is what is called the shear strain. The force divided by the area it is applied on is called the shear stress [4]. When a shear force is applied to a body the shear strain increases with time. When studying liquids their motion usually happens alongside a wall. In this case the velocity of the liquid changes depending on the distance from the wall. The variation of this velocity is what is called the shear rate. This and the shear stress needed to achieve a certain shear rate is how the viscosity is calculated.

If the viscosity is independent of the stress the fluid is said to be Newtonian. Blood is a fluid that has a shear thinning behaviour [5] which means that its viscosity decreases under shear strain.

Figure 2.1: Two shear forces in opposite directions caus-ing a shearcaus-ing deformation.

(8)

2.2 Elasticity

One important fundamental property of non-liquid materials is the Young’s modulus. Similar to the shear modulus the Young’s modulus is defined by the stress divided by the strain. It describes the elastic behaviour of a material or the modulus of elasticity. Using the dimensionless strain defined as the extension ratio (2.1) (the final length l divided by the original length L) and the force F per unit area A (2.2) (usually in N/m2 _{or Pascal), the Young’s modulus (E) can}

be defined as the ratio between them (2.3). Rubbers are an example of a material that have a relatively low Young’s modulus. Metals have a higher Young’s modulus [4].

= L l (2.1) σ = F A (2.2) E = σ (2.3)

2.3 Blood

In some regards blood behaves similarly to water, they both have a relatively low viscosity and therefore flow easily. But on a smaller scale there are more dynamics since it is not actually a pure fluid, rather it is a complex suspension of different particles in a fluid base [4]. These particles have different effects on the fluid behaviour of blood.

2.3.1 Plasma

Blood plasma is a non-Newtonian fluid [6] and is the component of blood that holds all the particles in suspension. Plasma consists for 90% of water and blood is made out of plasma for about 50% [4].

2.3.2 Red blood cells

Red blood cells are one of the main components in blood. The volume percentage, also known as the haematocrit, is usually between 30 and 50 percent depending on different factors such as sex and age [7]. Its main function is transferring oxygen from the lungs to body cells and carrying carbon dioxide back to the lungs. They typically live for around 120 days [8], have a mean diameter of 7.5µm [4], and a volume of 88 fL [2][9]. However, there is some distribution within these last two properties. For example when red blood cells age they lose volume [9].

Because red blood cells in blood all have a different age ranging from zero to 120 days, different RBC do not have the same dimensions. The variation of volume in red blood cells is also known as red blood cell distribution width (RDW). The RDW is calculated using the standard deviation of the mean cell volume. Typically the standard deviation is around 12 fL [2]. The RDW is then calculated as 12 fl/88 fl × 100 = 13.6%. The RDW for young adults ranges between 11.6 and 14.6 percent. It is important to note that this does not directly translate to the distribution in the diameter of the red blood cells since it is based on volume. A 13.6% increase in volume would roughly equate to a diameter increase of the cube root of 13.6%, 4.3%.

The shape of red blood cells is represented as a biconcave model in HemoCell and is shown in figure 2.2.

(9)

Figure 2.2: Red blood cell model as used in HemoCell [10].

2.3.3 White blood cells and platelets

Two other component of blood are white blood cells and platelets. White blood cells make up blood for about 0.7% of its volume [4], typically have a spherical shape, and vary in size with a diameter ranging from 7 to 30µm. Platelets make up 0.3% of the blood volume [4], are relatively small with a diameter of 2 to 3µm, and have a plate-like shape. These two cell types are also implemented in HemoCell but are not included in this research project.

(10)

(11)

CHAPTER 3

Implementation

3.1 Red blood cell diameter

The HemoCell framework was originally designed to receive one parameter for the different red blood cell properties. For example the diameter is set to be 7.82µm for all cells in the default configuration files. This is also the diameter that has been used in the model verification tests done by Z´avodsky et al. [1].

In-vivo red blood cells however, do not have one diameter value for all RBC. They actually vary in size [9]. To more accurately describe in-vivo blood flows in the HemoCell simulation model, a distribution can be applied on the RBC size using the RDW. To do this a random value needs to be picked from the distribution described by the RDW. This percentage can then be added to the original RBC volume.

Because the RBC volume distribution is a Gaussian distribution [9] these random values should not be generated using a simple uniform distribution. Some RBC will have a size that is bigger than the mean size plus RDW. The HemoCell framework already has a random number generator class that uses the drand48 function from the stdlib library which generates a random floating point number between 0.0 and 1.0. This class has to be extended with a new function that can generate numbers using a normal distribution.

The Box-Muller transform generates a pair of random numbers that are normally distributed using a pair of numbers between 0 and 1 that have a uniform distribution. If U1 and U2are two

random numbers between 0 and 1, two numbers with a normal distribution can be generated using 3.1 and 3.2 [11].

Z0= R cos(Θ) =p−2 ln U1cos(2πU2) (3.1)

Z1= R sin(Θ) =p−2 ln U1sin(2πU2) (3.2)

The numbers can then be multiplied by the preferred standard deviation, the RDW in this case, and added to the mean to create any normal distribution. To simulate this for a typical RBC size distribution a scaling factor can be used with the generated random numbers as its value. The mean would be 1.0 (a typical RBC) and the standard deviation is related to the RDW.

The scaling factor will scale the radius of the RBC. The RDW is a number based on the RBC volume and because volume and radius do not scale linearly the RBC radius should also be multiplied by a different value. Because the shape of RBC does not change much with size [3] a volume increase by a factor of 1.1 would mean that the diameter is actually only increased by a factor of √3_1.1.

(12)

3.2 Initial cell positions

Before the simulation can be started a list of cells is required. For each cell type being simulated a separate file is used. The files consist of six floating point numbers that describe the position in the simulation space and the three Euler angles yaw, pitch, and roll that tell how the cell is rotated. The first line of the file is a number that specifies exactly how many cells can be found in the file. The format of the file is described in listing 3.1. An example file containing three cells in random positions is shown in listing 3.2.

Listing 3.1: HemoCell position file format Ncells

x y z yaw p i t c h r o l l . . .

Listing 3.2: Example of HemoCell position file 3

9 9 . 4 8 7 3 0 . 7 1 6 4 0 7 3 2 . 9 7 2 8 6 9 . 8 5 9 2 5 2 . 9 9 4 3 8 4 . 9 8 5 6 9 . 8 9 0 2 9 1 5 . 3 9 1 6 2 6 . 1 7 2 8 −66.0462 4 8 . 5 2 2 7 1 6 7 . 5 9 6 2 1 . 2 6 7 8 4 9 . 7 6 1 2 1 1 . 2 2 4 1 2 5 . 7 0 7 2 7 . 5 8 7 3 1 4 8 . 3 4 8

These files are easily generated using the “packCells” tool. This tool is included with the Hemo-Cell package and initialises all cells at random positions and puts and simulates forces between these until they do not overlap. This packing process is not trivial and has been described before in another research paper [10]. During this process the cells are represented by ellipsoids that enclose the cell as closely as possible. This is done because when these ellipsoids do not overlap the cells can also not overlap in any way but ellipsoids are computationally simpler to check for collision. For a simulation with cells of varying size these ellipsoids need to be scaled individu-ally. To accommodate this a scaling factor needs to be added as a property for each cell and its ellipsoid. The factor is determined using a mean and standard deviation during the initialisation of every cell. A scaling factor of 1.0 means that the cell has the same size as the originally sized red blood cell. In the case of a distributed cell width of 13.6% this factor will be distributed randomly using a normal distribution with a standard deviation of √3

1.136 − 1 ≈ 0.0434 and a mean of 1.0. Using this distribution approximately 95% of the cells will have a scaling factor between 0.91 and 1.09 (µ − 2σ and µ + 2σ) as shown in figure 3.1.

The amount of cells that are being generated by the packCells tool is now also dependent on the mean size of the cells. In the original implementation a mean volume of 90 fl is assumed and the amount of cells generated is calculated by dividing the volume of the 3D space by the average cell size and multiplying it with the haematocrit. Because in some cases the mean will differ from 1.0, for example to simulate macrocytic or microcytic anemia, the new implementation needs to calculate the amount of generated cells by using the mean scaling factor. For a mean radius scaling factor of 0.95 the mean volume scaling factor is actually 0.953_{≈ 0.857. The mean}

volume for RBC in this case can be calculated as 90 fl × 0.953_{≈ 77.16 fl and the amount of cells}

generated for the same haematocrit value is increased by almost 17%.

A normal distribution of cell volume in red blood cells correlates with the physiological distribu-tion as described in [9]. This scaling factor is added to the posidistribu-tion file as a new cell property. The new format is described in listings 3.3 and 3.4.

(13)

Figure 3.1: Distribution of scaling factor.

1.96σ 0.91 1 1.09 scaling factor

probability

Listing 3.3: New HemoCell position file format Ncells

x y z yaw p i t c h r o l l s c a l i n g . . .

Listing 3.4: Example of new HemoCell position file 3

9 9 . 4 8 7 3 0 . 7 1 6 4 0 7 3 2 . 9 7 2 8 6 9 . 8 5 9 2 5 2 . 9 9 4 3 8 4 . 9 8 5 6 1 . 0 3 4 6 9 . 8 9 0 2 9 1 5 . 3 9 1 6 2 6 . 1 7 2 8 −66.0462 4 8 . 5 2 2 7 1 6 7 . 5 9 6 0 . 9 3 2 8 2 1 . 2 6 7 8 4 9 . 7 6 1 2 1 1 . 2 2 4 1 2 5 . 7 0 7 2 7 . 5 8 7 3 1 4 8 . 3 4 8 0 . 9 8 4

This format is read by HemoCell during the initialisation of the simulation. Some modifications to the HemoCell source code need to be made to include the scaling factor. One of the files that needs to be modified is readPositionsBloodCells.cpp which contains functions that are used to read the position file. The function getReadPositionsBloodCellsVector reads the position file and will have to output an extra vector for scaling next to the already existing vectors for cell-ids, angles, and positions. This vector is filled during the processing of the file by reading an extra floating point number for every line and pushing this number to the end of the vector.

This function is called during by another function processGenericBlocks which initialises all cells. In this function a mesh is created for all cells. Before this mesh can be used however it first needs to be scaled by the scaling factor that corresponds to the cell. This scaling factor is the one that was copied to a vector previously. Finally the cell gets rotated accordingly and put into simulation space by translating it to its position. It is important to first scale and then translate, otherwise the cell would not get the exact position but rather one that is factorised by the scaling factor.

3.3 Red blood cell rigidness

The rigidness of the red blood cells is defined in HemoCell using two properties, the link force coefficient (surface elasticity) and the bending force modulus which defines how easily an RBC

(14)

Table 3.1: Young’s modulus for different values of kLink.

Young’s modulus Surface elasticity (kLink, dimensionless)

69.3µN m−1 15 (default value in current HemoCell implementation) 43.7µN m−1 15.52

55.4µN m−1 19.86 67.1µN m−1 24.19

The elasticity for RBC in rabbits has been researched and an elasticity of around 55.4µN m−1 (Young’s modulus) with a standard deviation of 11.7µN m−1 has been found [9]. This research project also shows that cell age has no influence on the cell elasticity1_{, unlike the cell volume.}

For this reason the distributions for size and elasticity can be considered independent. By using the “materialTester” script included in the HemoCell toolset the Young’s modulus for different surface elasticity values can be numerically calculated. The surface elasticity of RBC in HemoCell is defined as a dimensionless variable called “kLink” with a default value of 15, which yields a Young’s modulus of 42.4µN m−1. To create a distribution for the surface elasticity that describes the distribution found in rabbits the dimensionless kLink values that correspond to the mean (55.4µN m−1) and the corresponding lower and upper bounds (55.4µN m−1− 11.7µN m−1 ₌

43.7µN m−1 and 55.4µN m−1+ 11.7µN m−1= 67.1µN m−1) need to be found.

By plotting the Young’s modulus with the corresponding HemoCell kLink value a linear distri-bution is found (figure 3.2). Applying linear regression on the numeric results gives the equations 1.787 + 2.7 × kLink = E and E−1.787_2.7 = kLink where E is Young’s modulus. Using this infor-mation the kLink values that correspond with the Young’s modulus found in the research paper can be calculated (table 3.1).

Because the relation between the dimensionless variable and elastic modulus is linear the normal distribution can be created using the kLink values found earlier. The mean will be 19.86 and the standard deviation can be calculated as 19.86 − 15.52 = 4.34. These two values can be used for kLink in HemoCell when running the simulation by distributing the variables using a normal distribution. Because RBC volume and elasticity are independent variables and the elasticity does not affect the initial positioning of the cells the value of kLink can be generated for each cell at runtime using a normal distribution as described in 3.1 instead of during the initialisation of the position file.

3.4 Simulation space

The simulation is run in a 3D space defined as a voxel space. Within this space a geometry can be defined which is usually done by loading an stl-file during the initialisation of the simulation. For modelling a blood vessel the geometry of a simple tube can be used as pictured in figure 3.3. HemoCell makes use of the Lattice Boltzmann methods for fluid simulation which operates on a grid. Because the stl-file describes the geometry surfaces by triangles it first needs to be voxelised before the simulation is run. HemoCell uses a Lattice Boltzmann library called Palabos that does this.

3.4.1 Actual haematocrit

Because during the initialisation of the cell positions using the packCell tool the geometry is not yet defined there will also be cells that either fall outside the geometry or overlap with the

1_{Another research however shows different degrees of transformation for young and old cells in shear flow [12].}

This might indicate that there actually is a relation between age and rigidness, which at the same time also means that there will be a non-random relation between diameter and rigidness. Treating these two variables as completely independent can be incorrect.

(15)

Figure 3.2: Relation between surface elasticity and elastic modulus. 0 10 20 30 40 50 0 30 60 90 120 150

surface elasticity (kLink ) [-]

Y oung’s mo dulus [µ N m − 1]

Value found numerically Fitted plot

Figure 3.3: 3D tube that is used to simulate the blood cells in.

geometry. A tube covers roughly 79% of its enclosing square prism2 _{which means that as the}

number of cells reach infinity 21% of all cells are expected to fall outside the geometry. However, cells that only overlap with the geometry edges are also marked for deletion. Since these cells are still partly inside the defined geometry the effective red blood cell volume inside the tube will go down when the cells overlapping with the geometry are removed.

Generating the initial cells first the number of required cells is calculated considering a mean volume of 90 fl. First the volume of the grid is divided by the mean volume of the red blood cells which gives the maximum amount of cells that could possibly fit in the defined space. This amount of cells is then multiplied by the haematocrit. In the case of a 100µm by 50 µm by 50µm grid and a haematocrit of 30%, 100 µm × 50 µm × 50 µm/90 fl × 30% = 833 red blood cells are generated. Around 21% of those are expected to fall outside the geometry and be marked for deletion. A simple test however shows that with the largest tube that fits inside the grid (a tube that has a length of 100µm and a radius of 25 µm) only 458 of the cells are not marked for deletion. This means that 45% of the cells are removed which is significantly more than the expected 21%.

Because more cells than just those outside the tube are deleted the haematocrit within the tube will decrease. The new haematocrit value can be found by calculating both the volume of the leftover red blood cells and the volume of the tube. The volume of the red blood cells

2_{This number can be reached by calculating the overlap of the biggest circle that fits inside a square and the}

square. For a square that is 10 by 10 a circle with a radius of 5 is the biggest circle that can fit inside. This gives a square with a surface area of 100 and a circle with a surface area of 25π. The overlap is 25π₁₀₀ × 100% ≈ 78.54. Since the tube and grid are simply prisms with a circle and square as a base the same will hold in the three

(16)

is expected to be around 458 × 90 fl = 41 220 fl for a mean RBC volume of 90 fl. The volume of the tube can be calculated as π252× 100 ≈ 196 349 fl. So the new haematocrit is actually

41220

196349× 100% ≈ 21%.

The same tube and a haematocrit of 40% yield 1111 generated cells of which only 625 are actually simulated which is still a removal of almost 44% making the haematocrit value inside the tube only 28.6%. This is something that needs to be taken into account when making a comparison to other models.

3.5 Velocity measurement

During the simulation csv-files containing data for each cell are generated at a configurable interval. These files contain the id of the cell, which is equal to the line number in the original position file minus one, and also its current position. This data can be used to calculate the average velocity within this interval.

For every interval HemoCell produces a csv file for each node the simulation runs on. So for a simulation that runs on 16 nodes 16 csv files are generated at every interval. First these files are merged into one file per interval. To analyse the remaining data a python program is written that reads every csv and stores the positions for each cell at every interval into memory. Because the cell positions that are written to the csv are in lattice units in the simulation grid they first need to be converted to micrometres. The lattice unit used for the simulations is 0.5µm so this is a simple division by two. After this the average speed over one interval can be calculated by calculating the distance within the euclidean space between the previous and current interval. Using this data the mean and standard deviation of the velocity over the whole or a part of the simulation can also be calculated.

By comparing this data to the original position file the cell size can be added to the cell velocities to see if any meaningful relationship can be found.

(17)

CHAPTER 4

Experiment

The original HemoCell model is compared to the new model with a heterogeneous red blood cell size distribution. The cells will be simulated in a 100µm long blood vessel with a radius of 25µm. The initial cell positions are generated using the haematocrit values 30% and 40% within a 100µm by 50 µm by 50 µm grid. The cells are simulated in the 3D tube pictured in figure 3.3 which represents a blood vessel. As described earlier in paragraph 3.4.1 the generation of cell positions with a 30% haematocrit value usually does not yield a haematocrit value within the tube that is 30%. This depends on how many cells completely fit inside the tube.

For the first three experiments the cells are simulated for 2 million iterations where each iteration is equal to 1 × 10−7 second. This makes a total simulation of 200 real-time milliseconds. The other three simulations are simulated for 120 real-time milliseconds. The simulation is run with a periodic boundary condition which means that when red blood cells get to the edge of the tube the cells are automatically moved back to the start of the tube keeping the same velocity and rotation. Cells that moves around 5000 micrometres per second will move through the whole tube around 6 times within the 1.2 million iterations. There are no cells removed or generated during the simulation.

This simulations are run with three different cell size distributions. One with all cells having a radius of 3.96µm, this is equivalent to the original HemoCell model. The second simulation will be a real-world situation where the radius is variant. The radius will be multiplied by a scaling factor found in in-vivo blood cell distributions (13.6%) with a mean of 1 and a standard deviation of √3

1.136 − 1 ≈ 0.0434. The third model will be an extrapolation of the real world model and uses an even bigger standard deviation of .129. An example of an initial state with a high varying size is shown in figure 4.1.

Table 4.1: Experiments run.

Experiment Initial haematocrit Cells generated Cells simulated Effective haematocrit Scaling standard deviation 1. 30% 833 458 20.99% 0 2. 30% 833 457 20.95% 0.0434 3. 30% 833 458 20.99% 0.1290 4. 40% 1111 622 28.51% 0 5. 40% 1111 621 28.46% 0.0434

(18)

(19)

CHAPTER 5

Results

5.1 Viscosity

During the simulation the bulk viscosity is calculated every 5000 iterations (half a millisecond in real-time). The positions of all individual cells are also written to a csv-file every 5000 iterations. The viscosity describes the force needed to move the flow. The viscosity for the first three experiments with a low haematocrit is shown in figure 5.1.

At the beginning of the simulation the viscosity rises steeply. This behaviour is exactly as described in previous research [1]. This is because at the beginning of the simulation the cells do not flow but deform. This process is shown in the figure 5.2. With the stabilisation of the cell deformation the viscosity also stabilises as visible from the graph.

During the initial starting up period which stops at around 33 milliseconds there seems to be a significant difference in how viscous the blood flows are. After this initial period however the graph plots start to intertwine. Around the 43rd millisecond of the simulation the three experiments show an almost identical viscosity property (1.32242, 1.32621, and 1.32719) which should be a good indicator that the flows have passed their initial phase. The average viscosity for the three flows at different times are shown in table 5.1. Although experiment 1 seems to have stabilised after 44 milliseconds, experiment 3 gets a sudden peak at around 30 milliseconds. This causes the average viscosity to be larger than that of experiment 2 in the beginning. Only looking at the first 20 milliseconds however, when all three experiments are still in the initial phase, the average relative viscosity of experiment 3 is smaller than that of experiment 2. Even though the average computed viscosity over the whole simulation is larger in experiment 1 than that of experiment 3 the viscosity is too unstable to say this is significant.

The surface area for every cell is also calculated and stored in the csv info file. The surface area can be used as a good indicator to see when the cells have fully deformed. The average surface area is shown figure 5.3 and confirms the stabilisation of the cell shearing effect to be around 44ms.

The bulk viscosity for the simulations with a larger haematocrit are shown in figure 5.4. During

Table 5.1: Calculated average relative viscosity.

Experiment Average relative viscosity after 42.5ms before 43ms before 20ms

1 1.3233 1.3404 1.3473

2 1.3221 1.3241 1.3275

(20)

Figure 5.1: Calculated bulk viscosity in the low-haematocrit simulations. 0 20 40 60 80 100 120 140 160 180 200 1.3 1.32 1.34 1.36 1.38 1.4 1.42 1.44 1.46 simulation time [ms] relativ e app aren t viscosit y [-]

Exp. 1, original, lower haematocrit Exp. 2, biological distribution, lower haematocrit

Exp. 3, extreme distribution, lower haematocrit

Figure 5.2: State of RBC at the initial state and after 3 and 30 milliseconds.

these simulations the initial period seems to be much longer but again the simulation with the original distribution still has the largest viscosity. After the initial period of around 100ms the computed viscosity starts to intertwine similarly to the low-haematocrit experiments.

5.2 Individual cell information

By analysing the csv files as described in paragraph 3.5 the velocity for individual cells can be measured. If the scaling property has a big effect on the cell velocity a trend should be visible in the scatter plots shown in figure 5.5. From the scatter plot there is no obvious trend. The bar diagrams in figure 5.6 do not show a significant influence on velocity by cell size either.

Something that does come to notice is the fact that the scatter plot does not show any cells with scaling factor of 1.3 or more. By analysing the histogram as shown figure 5.6 there does appear to be a normal distribution, but the mean is slightly shifted to the left. The mean scaling factor is actually calculated as 0.9468 and not the expected 1.0. A cause of this might be that bigger cells around the edge of the tube are more likely to clip with the geometry and get removed.

(21)

Figure 5.3: Average cell surface area during simulation. 0 10 20 30 40 50 60 70 80 90 100 126 128 130 132 134 simulation time [ms] a v erage cell v olume [µ m 2]

Average cell surface area (experiment 1) Average cell surface area (experiment 3)

Figure 5.4: Calculated bulk viscosity in the low-haematocrit simulations.

0 20 40 60 80 100 120 1.3 1.4 1.5 1.6 1.7 1.8 1.9 simulation time [ms] relativ e app are n t viscosit y [-]

Exp. 4, original, higher haematocrit Exp. 5, biological distribution, higher haematocrit

(22)

Figure 5.5: Average cell velocity over whole simulation, experiment 3. 0.6 0.8 1 1.2 1.4 0 2 4 6 8 scale factor [-] v elo cit y [µ m ms − 1] Velocity distribution 0.6 0.8 1 1.2 1.4 0 0.1 0.2 0.3 0.4 scale factor [-] standard deviation v elo cit y [µ m ms − 1]

Velocity standard deviation distribution

Figure 5.6: Average velocity for different scaling factors over whole simulation.

0.875-0.925 0.925-0.975 0.975-1.025 1.025-1.075 1.075-1.125 0 2 4 6 8

scaling factor, experiment 2

a v erage cell v elo cit y [µ m s − 1] 0.675-0.7250.725-0.7750.775-0.8250.825-0.8750.875-0.9250.925-0.9750.975-1.0251.025-1.0751.075-1.1251.125-1.1751.175-1.2251.225-1.275 0 2 4 6 8

scaling factor, experiment 3

a v erage cell v elo cit y [µ m s − 1]

(23)

Figure 5.7: Histogram showing frequency of cells at different sizes for experiment 3. 0.675-0.7250.725-0.7750.775-0.8250.825-0.8750.875-0.9250.925-0.9750.975-1.0251.025-1.0751.075-1.1251.125-1.1751.175-1.2251.225-1.275 0 20 40 60 80 scaling factor amoun t

cells within the tube is made. Figures 5.8 and 5.9 show the cell positions looking at only the y and z coordinates. Figure 5.8 shows the initial positioning at the start of the simulation which appears random which is expected when the random distribution of cells is functioning correctly. At 200ms the cell positions still appear random which shows that there is no obvious tendency for cells of a specific size to move outwards or to the centre.

The cells are also plotted as a function of distance from the centre and velocity in figure 5.10. As expected the velocity gets higher when the cells are closer to the centre. Using quadratic regression a plot is fitted. Calculating the average scaling factor above and below the fitted plot the values 0.952 (n = 205) and 0.942 (n = 243) were found suggesting that bigger cells might actually have a higher velocity. Similar results can be found with experiment 6. Experiment 2 and 4 however show almost identical average scaling factors for both below and above fitted plot. The slightly higher average speed could be because of bigger cells having a slightly higher chance to get removed at the edge than smaller cells, causing bigger cells to be more likely to be around the centre. The average distance for different scaling factors as shown in figure 5.11 does not show an obvious pattern either though, suggesting that the effect is quite small.

(24)

Figure 5.8: Cell positions in tube at 0ms, experiment 3. 0 10 20 30 40 50 60 70 80 90 100 0 20 40 60 80 100 y, lattice units [0.5µm] z, lattice units [0 .5 µ m] 0.7 0.8 0.9 1 1.1 1.2 1.3 scaling factor

Figure 5.9: Cell positions in tube at 200ms, experiment 3.

0 10 20 30 40 50 60 70 80 90 100 0 20 40 60 80 100 y, lattice units [0.5µm] z, lattice units [0 .5 µ m] 0.7 0.8 0.9 1 1.1 1.2 1.3 scaling factor

(25)

Figure 5.10: Average velocity at different distances from centre, experiment 3. 0 5 10 15 20 25 30 35 40 45 0 1 2 3 4 5 6 7 8

distance from centre, lattice units [0.5µm]

v elo cit y [µ m s − 1] 0.7 0.75 0.8 0.85 0.9 0.95 1 1.05 1.1 1.15 1.2 1.25 1.3 scaling factor

(26)

Figure 5.11: Average distance from centre for different cell scalings. 0.675-0.7250.725-0.7750.775-0.8250.825-0.8750.875-0.9250.925-0.9750.975-1.0251.025-1.0751.075-1.1251.125-1.1751.175-1.2251.225-1.275 30 32 34 36 38 40 scaling factor a v erage distance from cen tre, lattice units [0 .5 µ m]

(27)

CHAPTER 6

Conclusions / Discussion

This study gives insight on what the influence of a heterogeneous cell size distribution has on the viscous properties of blood flows and their red blood cells.

During the study the simulation run does not show a significant difference in the relative viscosity of the blood flow between the homogeneous and heterogeneous cell distribution during a long simulation. In the initial period however where the cells are still being effected by the shearing force there is a significant difference. Blood being a pulsatile flow this can be an interesting property. Further research could simulate a pulsatile flow in HemoCell to see what effect this exactly has on the individual cell and flow properties. A longer simulation can also be run to calculate the average viscosity with more significancy.

There also appears to be a slight tendency for bigger cells to have a higher velocity. This might have been partly caused by the fact that bigger cells are more likely to be removed from the edge of the tube. Further research could improve the packing algorithm such that geometry is already accounted for before the simulation starts. This would prevent cells from being removed and the haematocrit and mean cell size will be as expected instead of slightly lower.

There was a slightly incorrect distribution of cell sizes during the simulations which caused the mean cell size to be around 5% lower than expected. This in turn means that the actual simulated haematocrit is also 5% lower than expected for the experiment with a high variation. The experiments with a lower variation are less effected. This may have partly caused the initial viscosity peak to be higher for the blood flows with a low variation. A future study could rerun the simulation with a different scaling distribution. This can be done by using a higher haematocrit and bigger tube with more cells or by running the packing algorithm with a higher mean scaling factor value.

(28)

(29)

Bibliography

[1] G´abor Z´avodszky et al. “Cellular Level In-silico Modeling of Blood Rheology with An Improved Material Model for Red Blood Cells”. In: Frontiers in physiology 8 (2017), p. 563. [2] Neerja Vajpayee, Susan S Graham, and Sylva Bem. “Basic examination of blood and bone marrow”. In: Henry´s clinical diagnosis and management by laboratory methods 22 (2011), pp. 509–35.

[3] Michael K Obeng. “The American Heritage Stedman’s Medical Dictionary”. In: Journal of the National Medical Association 95.7 (2003), p. 634.

[4] P.R. Hoskins, P.V. Lawford, and B.J. Doyle. Cardiovascular Biomechanics. Springer Inter-national Publishing, 2017. isbn: 9783319464077.

[5] Pedram Tazraei, Alireza Riasi, and Behrouz Takabi. “The influence of the non-Newtonian properties of blood on blood-hammer through the posterior cerebral artery”. In: Mathe-matical biosciences 264 (2015), pp. 119–127.

[6] C Graf and J-P Barras. “Rheological properties of human blood plasma - a comparison of measurements with three different viscometers”. In: Experientia 35.2 (1979), pp. 224–225. [7] CK Cheng et al. “Complete blood count reference interval diagrams derived from NHANES III: stratification by age, sex, and race.” In: Laboratory hematology: official publication of the International Society for Laboratory Hematology 10.1 (2004), pp. 42–53.

[8] Christopher J Lindsell et al. “A method for the continuous calculation of the age of labeled red blood cells”. In: American journal of hematology 83.6 (2008), pp. 454–457.

[9] Richard E Waugh et al. “Rheologic properties of senescent erythrocytes: loss of surface area and volume with red blood cell age”. In: Blood 79.5 (1992), pp. 1351–1358.

[10] G´abor Zavodszky et al. “Hemocell: a high-performance microscopic cellular library”. In: Procedia Computer Science 108 (2017), pp. 159–165.

[11] George EP Box, Mervin E Muller, et al. “A note on the generation of random normal deviates”. In: The annals of mathematical statistics 29.2 (1958), pp. 610–611.

[12] Conrad Pfafferott, Gerard B Nash, and Herbert J Meiselman. “Red blood cell deformation in shear flow. Effects of internal and external phase viscosity and of in vivo aging”. In: Biophysical journal 47.5 (1985), pp. 695–704.

The Effect of a Non-Homogenous Erythrocytes Properties Distribution in In-Silico Blood Flows

Bachelor Informatica