A method for simulating the process of diabetic retinopathy

Abstract

The present invention discloses a method for simulating the process of diabetic retinopathy, including constructing a two-dimensional discrete coordinate system and defining a vascular index function that describes the cell grid properties of the cell grid; updating the blood glucose concentration, oxygen concentration, fibronectin concentration, and growth factor concentration of each cell grid at preset intervals t; predicting the movement direction of each blood vessel tip within the current time period t based on the fibronectin concentration and growth factor concentration; calculating a first predicted extension distance and a second predicted extension distance based on the movement direction within the current time period t and the movement direction within the previous time period t; and updating the cell grid properties of each cell grid in the two-dimensional discrete coordinate system based on the movement direction within the current time period t, the vessel tip's existence time, and the calculated movement extension distance of the vessel tip. The present invention incorporates oxygen concentration and blood glucose concentration into a model for simulating blood vessel extension, achieving accurate simulation prediction of diabetic retinopathy.

Description

Technical Field

[0001] The present invention relates to the technical field of biological simulation, and in particular to a method for simulating the process of diabetic retinopathy. Background Art

[0002] Diabetic retinopathy (DR) is a common complication in diabetic patients. The retina is the most light-sensitive tissue in the eye, and studies have shown that its health depends on a stable blood and oxygen supply. However, long-term hyperglycemia can damage retinal blood vessels, causing retinal microcirculation disorders and vascular lesions, leading to retinal hypoxia, hemorrhage, exudation, and other problems, destroying the retina's stable blood and oxygen supply environment and causing retinopathy. It can be seen that the high blood sugar environment of diabetic patients and the oxygen concentration in the high blood sugar environment are important pathological factors leading to diabetic retinopathy. Establishing a mathematical model that includes the influence of oxygen concentration and blood sugar is of great significance for a deeper understanding of the mechanism of diabetic retinopathy and the development of new treatments.

[0003] The simulation of diabetic retinopathy has important scientific significance and clinical application value. The realization of simulation of diabetic retinopathy process helps to deeply understand the disease mechanism, optimize treatment plans, discover new treatment methods, and provide support for personalized medicine and medical education.

[0004] Currently, research has been devoted to understanding the pathological mechanisms of diabetic retinopathy and developing effective treatments. Among them, the role of growth factors (such as VEGF) in retinal angiogenesis has been widely studied. However, for diabetic patients, a high blood sugar environment is a major symptom. In the process of diabetic retinopathy research, the prediction of the retinal angiogenesis process without considering the blood sugar level of diabetic patients is far from the actual situation. Summary of the Invention

[0005] In order to solve the problem that existing research efforts to understand the pathological mechanism of diabetic retinopathy and develop effective treatment methods do not consider the blood glucose level of diabetic patients, resulting in a significant difference between the predicted results and the actual situation, an embodiment of the present invention provides a method for simulating the diabetic retinopathy process.

[0006] The specific technical solutions are as follows:

[0007] An embodiment of the present invention provides a method for simulating diabetic retinopathy, comprising the following steps:

[0008] Step S100: constructing a two-dimensional discrete coordinate system and defining a blood vessel index function n(x, y, t) that describes and records the cell grid attributes of each cell grid in the two-dimensional discrete coordinate system at each preset time t; wherein the cell grid attributes include the blood vessel, the blood vessel tip, and the environment;

[0009] Step S200: Initializing the blood glucose concentration and oxygen concentration of each unit grid in the two-dimensional discrete coordinate system according to the patient's examination results, and initializing the unit grid attributes, fibronectin concentration, and growth factor concentration of each unit grid in the two-dimensional discrete coordinate system according to a preset initial state;

[0010] Step S300: updating the blood glucose concentration, oxygen concentration, fibronectin concentration, and growth factor concentration of each unit grid at a preset time interval t based on the unit grid's current blood glucose concentration, the governing equation for the blood glucose concentration, the current oxygen concentration, the governing equation for the oxygen concentration, the current unit grid attributes, the current fibronectin concentration, the governing equation for the fibronectin concentration, the current growth factor concentration, and the growth factor diffusion reaction equation;

[0011] Step S400, predicting the movement direction of each blood vessel tip in the two-dimensional discrete coordinate system within the current time period t based on the fibronectin concentration of the unit grid after time t, the growth factor concentration after time t, and the blood vessel tip movement direction control equation;

[0012] Step S500: Predicting a first motion velocity of each blood vessel tip in a two-dimensional discrete coordinate system when moving in a first motion direction based on a governing equation for the fibronectin concentration at the blood vessel tip after time t, the first motion direction, and the motion velocity of the blood vessel tip; and predicting a second motion velocity of each blood vessel tip in a two-dimensional discrete coordinate system when moving in a second motion direction based on a governing equation for the fibronectin concentration at the blood vessel tip after time t, the second motion direction, and the motion velocity of the blood vessel tip; wherein the first motion direction is the motion direction of the blood vessel tip during the current time period t; and the second motion direction is the motion direction of the blood vessel tip during the previous time period t.

[0013] Step S600: Calculate a first predicted extension distance of each blood vessel tip in a first motion direction in a two-dimensional discrete coordinate system based on the discrete coordinates (x0, y0) of the blood vessel tip and a first motion velocity, and calculate a second predicted extension distance of each blood vessel tip in a second motion direction in the two-dimensional discrete coordinate system based on the discrete coordinates (x0, y0) of the blood vessel tip and a second motion velocity; wherein the first predicted extension distance or the second predicted extension distance is the total number of unit grids through which the blood vessel tip moves when it starts at the discrete coordinates (x0, y0) of the unit grid to which the blood vessel tip belongs in the first motion direction or the second motion direction at the first motion velocity or the second motion velocity and reaches the end point coordinates (x1, y1) after a preset time t.

[0014] Step S700: Perform an extension simulation on each blood vessel tip in the two-dimensional discrete coordinate system within the current time period t based on the first movement direction of the blood vessel tip, the existence time of the blood vessel tip, the second movement direction of the blood vessel tip, the first movement speed, the second movement speed, the first predicted extension distance, and the second predicted extension distance. The unit grid attributes and the existence time of each unit grid in the two-dimensional discrete coordinate system are updated as the extension simulation of each blood vessel tip is performed. The simulation time is updated when the extension simulation of all blood vessel tips in the two-dimensional discrete coordinate system is completed. The extension simulation includes blood vessel extension, blood vessel branching that occurs during blood vessel extension, and blood vessel fusion that occurs during blood vessel extension.

[0015] Optionally, step S300 includes the following steps:

[0016] Predict the blood glucose concentration of each unit grid after time t based on the current blood glucose concentration of the unit grid and the control equation of the blood glucose concentration; and save the information under the corresponding unit grid;

[0017] Predict the oxygen concentration of each unit grid after time t based on the current oxygen concentration of the retinal unit grid, the blood glucose concentration after time t, and the control equation of the oxygen concentration; and save it under the information of the corresponding unit grid;

[0018] Predicting the fibronectin concentration of each retinal unit grid after time t based on the current unit grid properties of the retinal unit grid, the current fibronectin concentration and the control equation of the fibronectin concentration, and saving it under the information of the corresponding unit grid;

[0019] According to the current oxygen concentration, growth factor concentration and growth factor diffusion reaction equation of the unit grid, the growth factor concentration of each unit grid after time t is predicted, and the information of the corresponding unit grid is saved.

[0020] Optionally, the control equation of the oxygen concentration is as follows:

[0021]

[0022] in:

[0023] c o is the current oxygen concentration of the retinal cell grid;

[0024] D o is the diffusion coefficient of oxygen in retinal blood vessels, which is a fixed constant;

[0025] F(c o ) is the oxygen consumption rate equation of the retinal unit grid:

[0026]

[0027] A is the normal oxygen consumption rate of the unit grid, which is a fixed constant;

[0028] B is the constant of oxygen consumption of the unit grid under the influence of blood glucose;

[0029] c olim The lowest oxygen concentration for normal physiological activities of the unit grid is a fixed constant. Unit grids below this constant value are in an oxygen-deficient state.

[0030] G(ρ) is the amount of oxygen transported from the surrounding blood vessels in the unit grid. Vascular transport is the source of oxygen for the retinal unit grid. Based on the current vascular density, the calculation equation for G(ρ) is as follows:

[0031]

[0032] α o is the oxygen transport rate constant in blood vessels;

[0033] C is the transport efficiency constant under the influence of blood glucose;

[0034] ρ(x, y, t0) is the blood vessel density of the retinal unit grid (x, y) at the current time t0, specifically the number of unit grids with blood vessels in the current unit grid (x, y) and the nine unit grids adjacent to the current unit grid.

[0035] Optionally, the control equation for the fibronectin concentration is:

[0036]

[0037] in:

[0038] c f is the current fibronectin concentration of the retinal unit grid; the initial value has been set in step S200;

[0039] ω is the production rate of fibronectin in endothelial cells, which is a fixed constant;

[0040] is the uptake rate of fibronectin by endothelial cells, which is a fixed constant.

[0041] Optionally, the growth factor diffusion reaction equation is as follows:

[0042]

[0043] in:

[0044] c υ is the current growth factor VEGF concentration of the retinal unit grid; the initial value is the steady-state value under normal conditions;

[0045] D is the diffusion coefficient matrix of growth factor VEGF, specifically:

[0046]

[0047] d is a fixed value;

[0048] λ * is the half-life degradation rate, which is a constant value;

[0049] λ is the consumption rate of growth factor VEGF by endothelial cells, which is a constant value;

[0050] η is the production efficiency of growth factor VEGF under the control of oxygen concentration, which is a constant value;

[0051] S(C o ) is the control equation for generating the growth factor concentration under the influence of the current oxygen concentration of the retinal unit grid; if the oxygen concentration of the unit grid is C lim The above does not produce growth factor VEGF due to lack of oxygen. If the oxygen concentration of the unit grid is less than C lim VEGF will be produced due to lack of oxygen, so the equation is as follows:

[0052]

[0053] Optionally, the control equation for the movement direction of the blood vessel tip is as follows:

[0054]

[0055] in:

[0056] is the chemotactic parameter of growth factor VEGF;

[0057] is the fibrin chemotaxis parameter;

[0058] K is the extracellular matrix anisotropy tensor matrix, which is an intrinsic property of the matrix space of each cell grid (x, y)∈Ω and is independent of time. It is non-uniform and anisotropic in space. This property means that the extension direction of the tip cell is not necessarily consistent with the direction of the growth factor concentration gradient or the fibronectin concentration gradient; the governing equation is:

[0059]

[0060] k a is the anisotropy parameter, representing the resistance to movement along a certain direction;

[0061] k cond Heterogeneity parameters;

[0062] υ x and υ y is the component of the unit vector υ that describes the orientation of the extracellular matrix at the unit grid (x, y). At each unit grid (x, y), is a regular equation, randomly generating υ x and υ y , and fixed at all times t; for each unit grid (x, y) ∈ Ω, k cond and k a Also randomly generated.

[0063] Optionally, the control equation for the velocity of the blood vessel tip is as follows:

[0064]

[0065] in:

[0066] P is the velocity of the blood vessel tip;

[0067] k p is the rate of vascular tip division, which is related to the concentration of growth factor VEGF c υ Related.

[0068]

[0069] l e is a constant, the length of the endothelial cell;

[0070] t c The cell division cycle is determined by the concentration of growth factor VEGF c υ Determine that the control equation is as follows:

[0071]

[0072] τ is a constant, the cell proliferation time coefficient;

[0073] is the average growth factor (VEGF) concentration; the sum of the growth factor concentrations of all unit grids in the two-dimensional discrete coordinate system at the current time is divided by the total number of unit grids, and is calculated once for each time period t;

[0074] c υ is the current concentration of growth factor VEGF;

[0075] c υlim is the threshold concentration of growth factors.

[0076] Optionally, step S600 includes the following steps:

[0077] According to the velocity of the blood vessel tip in the x and y directions (x p ,y p ) and the discrete coordinates (x0, y0) of the blood vessel tip to calculate the end point coordinates (x1, y1) of the cell grid this extension; where x1 = x0 + x p , y1=y0+y p ; The movement speed is the first movement speed or the second movement speed;

[0078] The Bresenham algorithm is used to calculate the discrete coordinates of the unit grids passed on the straight line connecting the discrete coordinates (x0, y0) of the blood vessel tip to the end point coordinates (x1, y1) of the current extension of the blood vessel tip, and the discrete coordinates of all the passed unit grids are used as the predicted extension distance of the blood vessel tip (x0, y0) within the current time period t; the predicted extension distance is the first predicted extension distance or the second predicted extension distance.

[0079] Optionally, the conditions for the blood vessel branching when the blood vessel tip extends include:

[0080] The vessel existence time of the vessel tip is greater than a preset vessel tip existence time threshold τ a ;

[0081] The unit velocity vector of the blood vessel tip during this time period t The angle change between the angle between the positive direction of the x-axis of the coordinate system and the angle between the unit velocity vector of the blood vessel tip and the positive direction of the x-axis of the coordinate system in the last time period t is greater than the preset angle threshold τ u ;

[0082] When the two conditions for vascular branching are met at the same time, the process of vascular branching at the vascular tip is simulated as follows:

[0083] Starting from the tip of the blood vessel where the blood vessel branch occurs, the blood vessel is extended in the second movement direction according to the second predicted extension distance and in the first movement direction according to the first predicted extension distance;

[0084] The conditions for vascular fusion to occur when the vascular tip extends include:

[0085] Among the cell meshes included in the blood vessel tip prediction extension distance, there are cell meshes whose cell mesh attributes are blood vessels and / or blood vessel tips.

[0086] When the vascular fusion conditions are met, the vascular fusion process at the vascular tip is simulated including:

[0087] According to the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point on the corresponding predicted extension distance, the blood vessel tip is sequentially fused with the blood vessels and blood vessel tips on the predicted extension distance, and the blood vessel tip is extended according to the actual extension distance after the blood vessel fusion, until the extension simulation of the blood vessel tip in this time period t is completed;

[0088] The blood vessel tip extending the blood vessel comprises:

[0089] Each blood vessel tip in the two-dimensional discrete coordinate system is sequentially extended within the current time period t by simulation, wherein the blood vessel tip extension simulation within the current time period t is:

[0090] Determine whether vascular branching occurs at the vascular tip;

[0091] If no branching occurs at the vascular tip, it is determined whether vascular fusion occurs at the vascular tip within the first predicted extension distance. If no vascular fusion occurs, the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids within the first predicted extension distance. The current vascular tip existence time of the vascular tip is used as the initial value of the timer for the last unit grid within the actual extension distance. The vascular tip existence time of the vascular tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to vascular tip, and the timer is started to record the vascular tip existence time from the initial value. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessels, and the extension simulation of the vascular tip within the current time period t is terminated. If vascular fusion occurs, the vascular tip is sequentially fused with the blood vessels and the vascular tip within the first predicted extension distance in the order from the discrete coordinates of the vascular tip to the end point coordinates within the first predicted extension distance. The vascular tip is extended according to the actual extension distance of the vascular tip after vascular fusion until the extension simulation of the vascular tip within the current time period t is terminated.

[0092] If a blood vessel branch occurs at the tip of the blood vessel, the tip of the blood vessel is extended in the second movement direction according to the second predicted extension distance and the tip of the blood vessel is extended in the first movement direction according to the first predicted extension distance; wherein:

[0093] The extending of the blood vessel tip in the first movement direction according to the first predicted extension distance specifically includes:

[0094] determining whether vascular fusion occurs at the tip of the blood vessel at the first predicted extension distance;

[0095] If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids in the first predicted extension distance, the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped, the unit grid attribute of the last unit grid in the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero, and the unit grid attributes of the remaining unit grids in the actual extension distance are changed to blood vessel, and the extension simulation of the blood vessel tip in this time period t is ended;

[0096] If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips at the first predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point at the first predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed;

[0097] The extending of the second movement direction of the blood vessel tip according to the second predicted extension distance specifically includes:

[0098] determining whether vascular fusion occurs at the tip of the blood vessel at the second predicted extension distance;

[0099] If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids on the second predicted extension distance, and the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid on the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero. The unit grid attributes of the remaining unit grids on the actual extension distance are changed to blood vessels, and the extension simulation of the blood vessel tip in this time period t is ended.

[0100] If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips on the second predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point on the second predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed.

[0101] Optionally, fusing the blood vessel tip with the blood vessel and the blood vessel tip at the predicted extension distance, and extending the blood vessel tip according to the actual extension distance of the blood vessel tip after the fusion includes:

[0102] If the vessel fused with the blood vessel tip is within the predicted extension distance, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids before the unit grid to which the blood vessel tip belongs within the predicted extension distance, the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped, and the unit grid attributes of all unit grids within the actual extension distance are changed to blood vessels, thereby ending the extension simulation of the blood vessel tip within the current time period t.

[0103] If the blood vessel tip that undergoes vascular fusion is the blood vessel tip at the predicted extension distance, the blood vessel fusion is performed according to the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the predicted extension distance, specifically:

[0104] If the angle between the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the corresponding predicted extension distance is greater than or equal to 90 degrees, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids before the unit grid to which the blood vessel fusion occurs at the predicted extension distance. The blood vessel tip existence time of the blood vessel tip and the blood vessel tip at the predicted extension distance is reset to zero and the timing is stopped. The unit grid attributes of all unit grids at the actual extension distance are changed to blood vessels, and the extension simulation of the blood vessel tip in this time period t is terminated.

[0105] If the angle between the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the corresponding predicted extension distance is less than 90 degrees, if the blood vessel tip is randomly selected for extension, the blood vessel tip at the predicted extension distance stops extending; when there is no blood vessel or blood vessel tip after the blood vessel tip at the predicted extension distance that merges with the blood vessel tip at the corresponding predicted extension distance, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids included in the predicted extension distance of the blood vessel tip, and the blood vessel tip existence time of the blood vessel tip and the blood vessel tip at the predicted extension distance that merges with it is cleared and the timing is stopped, the unit grid attribute of the last unit grid at the actual extension distance is changed to the blood vessel tip and the timer is started to count from the initial value Recording the existence time of the blood vessel tip and changing the unit grid attributes of the remaining unit grids in the actual extension distance to blood vessels, ending the extension simulation of the blood vessel tip in this time period t; when there are still blood vessels and / or blood vessel tips after the blood vessel tip in the predicted extension distance that merges with the blood vessel tip in the corresponding predicted extension distance, the blood vessel tip continues to fuse with the blood vessel and blood vessel tip in the predicted extension distance, and extends according to the actual extension distance of the blood vessel tip after blood vessel fusion, until the extension simulation of the blood vessel tip in this time period t is ended; wherein, the initial value of the timer of the last unit grid in the actual extension distance is the existence time of the blood vessel tip before the blood vessel tip is cleared if no blood vessel branch occurs, and the initial value is zero if a blood vessel branch occurs;

[0106] If the vascular tip is randomly selected to stop extending, the vascular tip is extended at the predicted extension distance where the vascular tip merges with the vascular tip; the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids before the unit grid to which the vascular tip belongs at the predicted extension distance where the vascular tip merges with the vascular tip at the corresponding predicted extension distance, and the vascular tip existence time of the vascular tip is cleared and the timing is stopped, and the attributes of all unit grids at the actual extension distance are changed to blood vessels, ending the extension simulation of the vascular tip within this time period t.

[0107] The embodiment of the present invention provides a method for simulating the process of diabetic retinopathy, which has the following technical effects:

[0108] (1) The impact of blood glucose levels in diabetic patients on diabetic retinopathy is fully considered, and oxygen concentration and blood glucose concentration are incorporated into the simulated vascular extension model. This not only achieves accurate simulation predictions, but also allows for a more comprehensive simulation of the occurrence and development of diabetic retinopathy. Understanding how these factors interact to cause retinal vascular abnormalities helps reveal the root cause of the disease;

[0109] (2) Current treatments for diabetic retinopathy include laser therapy and anti-VEGF injections. By dynamically simulating the diabetic retinopathy process at different stages, this application can predict disease progression and the effects of different treatments, optimize treatment strategies at different stages, and improve treatment outcomes.

[0110] (3) Considering the blood sugar control level and oxygen supply status of individual patients, a personalized simulation model is established for different patients to achieve personalized disease prediction and treatment plan optimization;

[0111] (3) By simulating the regulatory effects of different drugs on oxygen concentration and blood sugar levels, potential new drugs can be screened and their efficacy can be evaluated, accelerating the new drug development process.

[0112] Other features and advantages of the present invention will be described in the following description, and in part will become apparent from the description, or will be understood by practicing the present invention. The purpose and other advantages of the present invention can be realized and obtained by the structures particularly pointed out in the written description and the accompanying drawings.

[0113] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0114] The accompanying drawings are used to provide a further understanding of the present invention and constitute a part of the specification. Together with the embodiments of the present invention, they are used to explain the present invention and do not constitute a limitation of the present invention. In the accompanying drawings:

[0115] Figure 1 A schematic flow chart of a method for simulating the process of diabetic retinopathy provided by an embodiment of the present invention;

[0116] Figure 2 Flowchart of step S700 in an embodiment of the present invention. DETAILED DESCRIPTION

[0117] The preferred embodiments of the present invention are described below with reference to the accompanying drawings. It should be understood that the preferred embodiments described herein are only used to illustrate and explain the present invention, and are not used to limit the present invention.

[0118] Oxygen concentration: Studies have shown that retinal hypoxia is a key factor in causing abnormal angiogenesis. However, the specific distribution of oxygen concentration and its dynamic changes in the mechanism of action in diabetic retinopathy remain unclear.

[0119] Blood glucose levels: Hyperglycemia is a key characteristic of diabetes. It is known to cause endothelial cell dysfunction and increase oxidative stress. However, the specific mechanisms by which blood glucose levels exacerbate retinopathy by affecting angiogenesis and oxygen supply require further investigation.

[0120] Comprehensive Model Construction: This study will, for the first time, incorporate oxygen concentration and blood glucose levels into the mathematical model of angiogenesis, constructing a more comprehensive model of diabetic retinopathy. By simulating the retinal angiogenesis process under different oxygen and blood glucose conditions, the mechanisms of their interaction will be revealed.

[0121] In diabetic patients, hyperglycemia is the main characteristic of diabetes. Hyperglycemia can lead to damage to retinal microvessels and aggravate hypoxia, thereby triggering the following series of biological reactions:

[0122] High blood sugar triggers inflammation and cell damage: Long-term high blood sugar can induce abnormal metabolic products in cells, triggering an inflammatory response, leading to retinal tissue cell damage, blood vessel blockage and constriction, and affecting blood supply;

[0123] Exacerbated hypoxia: As hyperglycemia affects the normal function of microvessels, retinal tissue hypoxia is further aggravated;

[0124] Release of angiogenic factors: Hypoxia stimulates retinal tissue to release angiogenic factors, such as vascular endothelial growth factor (VEGF), which promotes the growth of new blood vessels;

[0125] Increased angiogenesis: Under the action of angiogenic factors, the retinal tissue begins to abnormally proliferate new blood vessels. These new blood vessels are unstable and prone to rupture and bleeding, posing a serious threat to vision.

[0126] In summary, in order to simulate the complex phenomena and system behaviors in the process of diabetic retinopathy, help researchers understand and predict the complexity of the diabetic retinopathy process, analyze how blood glucose levels aggravate retinopathy by affecting angiogenesis and oxygen supply, and the role of oxygen concentration distribution and its dynamic changes in diabetic retinopathy, this application incorporates oxygen concentration and blood glucose concentration into the mathematical model of angiogenesis for the first time to construct a more comprehensive diabetic retinopathy model. By simulating the retinal angiogenesis process under different oxygen and blood glucose conditions, the interaction mechanism is revealed. Figure 1 As shown, the embodiment of the present application provides a method for simulating the process of diabetic retinopathy, which includes the following steps:

[0127] Step S100: construct a two-dimensional discrete coordinate system (x, y), and define a vascular index function n(x, y, t) that describes and records the attributes of each retinal unit grid in the two-dimensional discrete coordinate system at each preset time t; the unit grid attributes include blood vessels, blood vessel tips, and environment; it can be seen that the vascular index function n(x, y, t) in the embodiment of the present application is actually a binary value function with a value of 0, 1, or 2, indicating the attribute of the unit grid (x, y) at time t. For example, n(x, y, t) is 0, indicating that the unit grid attribute at time t is the environment, and n(x, y, t) is 1, indicating that the unit grid attribute at time t is the environment. When n(x, y, t) is 2, it means that the cell grid attribute at time t is the blood vessel tip. A binary value function is used as the vascular network morphology descriptor, which captures grid information more accurately than the traditional endothelial cell density function. The blood vessel tip is a state before it becomes a blood vessel. The biological processes such as blood vessel tip extension, branching and fusion are used to simulate the generation process of retinal blood vessels. The cell grid attribute is environment, which is a space for simulating blood vessel growth. The cell grid with the cell grid attribute of environment may become a cell grid with the cell grid attribute of blood vessel or blood vessel tip as the blood vessel tip extends, branches and fuses.

[0128] Step S200: Start simulation for a specific patient and initialize the properties of each cell grid in the two-dimensional discrete coordinate system, including blood glucose concentration, oxygen concentration, fibronectin concentration, and growth factor concentration. To ensure that the simulation results are closer to the patient's actual condition and closer to actual development, the blood glucose and oxygen concentrations are set based on the patient's examination results. The initial values are directly obtained from the patient's most recent test results. The initial values of fibronectin and growth factor concentrations are the steady-state values under normal conditions. The initial values of blood glucose concentration, oxygen concentration, and fibronectin concentration are the same for each cell grid. Over time, the blood glucose concentration, oxygen concentration, and fibronectin concentration of each cell grid change individually, and the blood glucose concentration, oxygen concentration, and fibronectin concentration vary from cell grid to cell grid.

[0129] In an embodiment of the present application, when initializing the attributes of each unit grid in a two-dimensional discrete coordinate system, the unit grid attribute of at least one unit grid is the blood vessel tip, and the attributes of the remaining unit grids are the environment or blood vessels; the position of the unit grid whose unit grid attribute is the blood vessel tip is not restricted here, but in order to better demonstrate the effect of blood glucose concentration on angiogenesis over time, when initializing in this embodiment, the number of unit grids whose unit grid attribute is the blood vessel tip is much smaller than the number of unit grids whose unit grid attribute is the environment, and the number of unit grids whose unit grid attribute is the blood vessel can be zero or less than the number of unit grids whose unit grid attribute is the environment.

[0130] Step S300: Update the blood glucose concentration, oxygen concentration, fibronectin concentration, and growth factor concentration of each unit grid at preset time intervals t based on the current blood glucose concentration of the unit grid, the control equation of the blood glucose concentration, the current oxygen concentration, the control equation of the oxygen concentration, the current unit grid properties, the current fibronectin concentration, the control equation of the fibronectin concentration, the current growth factor concentration, and the growth factor diffusion reaction equation, and save them under the information of the corresponding unit grid.

[0131] In the embodiment of the present application, step S300 includes the following steps:

[0132] The first step is to predict the blood glucose concentration of each unit grid after time t based on the current blood glucose concentration of the retinal unit grid and the control equation of the blood glucose concentration; and save it under the information of the corresponding unit grid.

[0133] The control equation for the retinal blood glucose concentration in the embodiment of the present application is input according to the patient's specific conditions, and is as follows:

[0134]

[0135] in:

[0136] Cg is the current blood glucose concentration of the retinal unit grid; the initial value has been set in step S200;

[0137] C g ( t ) is the current blood glucose concentration Cg based on the retinal unit grid, and the blood glucose concentration after time t of the retinal unit grid;

[0138] F ( t ) is the increase in blood glucose concentration at time t after food intake;

[0139] I ( t ) is the increase in insulin concentration after time t, and the initial value of insulin concentration is obtained from the patient's test;

[0140] H ( t ) is the increase in glucagon concentration after time t, and the initial value is obtained from the patient's test;

[0141] k1 is the influence coefficient of insulin concentration on blood glucose concentration, which is a fixed constant;

[0142] k2 is the coefficient of influence of glucagon concentration on blood glucose concentration, which is a fixed constant;

[0143] L ( t )is the blood glucose concentration consumed by the liver effect after time t.

[0144] Step 2: Predict the oxygen concentration of each unit grid after time t based on the current oxygen concentration of the retinal unit grid, the blood glucose concentration after time t, and the control equation of the oxygen concentration; and save it under the information of the corresponding unit grid;

[0145] In the embodiment of the present application, the change in oxygen concentration over time is mainly due to the generation, diffusion, and consumption of oxygen. In the retina, oxygen mainly comes from the retinal blood vessels. The simulation needs to consider the process of oxygen diffusion from blood vessels to retinal tissue. Therefore, the control equation of oxygen concentration in the retina is as follows:

[0146]

[0147] in:

[0148] c o is the current oxygen concentration of the retinal unit grid, the initial value of which has been set in step S200;

[0149] D o is the diffusion coefficient of oxygen in retinal blood vessels, which is a fixed constant;

[0150] F(c o ) is the oxygen consumption rate equation of the retinal unit grid:

[0151]

[0152] A is the normal oxygen consumption rate of the unit grid, which is a fixed constant;

[0153] B is the constant of oxygen consumption of the unit grid under the influence of blood glucose;

[0154] c olim The lowest oxygen concentration for normal physiological activities of the unit grid is a fixed constant. Unit grids below this constant value are in an oxygen-deficient state.

[0155] G(ρ) is the amount of oxygen transported from the surrounding blood vessels in the unit grid. Vascular transport is the source of oxygen for the retinal unit grid. Based on the current vascular density, the calculation equation for G(ρ) is as follows:

[0156]

[0157] α o is the oxygen transport rate constant in blood vessels;

[0158] C is the transport efficiency constant under the influence of blood glucose;

[0159] ρ(x, y, t0) is the vascular density of the retinal unit grid (x, y) at the current time t0, specifically the number of unit grids with blood vessels in the current unit grid (x, y) and the nine unit grids adjacent to the current unit grid;

[0160] Step 3: predict the fibronectin concentration of each retinal unit grid after time t based on the current unit grid properties of the retinal unit grid, the current fibronectin concentration and the control equation of the fibronectin concentration, and save it under the information of the corresponding unit grid;

[0161] Fibronectin is attached to the extracellular matrix and is an important substance for angiogenesis. The increase in fibronectin concentration depends on the presence of endothelial cells, while the decrease in fibronectin concentration depends on the uptake of fibronectin by endothelial cells. It can be seen that whether it increases or decreases depends on the unit grid attribute being a blood vessel, that is, it is related to the value of the vascular index function n(x, y, t), that is, the unit grid attribute of the unit grid (x, y) at t=t0. Therefore, the control equation of the fibronectin concentration in the embodiment of the present application, which does not include the diffusion term, is:

[0162]

[0163] in:

[0164] c f is the current fibronectin concentration of the retinal unit grid; the initial value has been set in step S200;

[0165] ω is the production rate of fibronectin in endothelial cells, which is a fixed constant;

[0166] μ is the uptake rate of fibronectin by endothelial cells and is a fixed constant.

[0167] Step 4: predict the growth factor concentration of each unit grid after time t based on the current oxygen concentration, growth factor VEGF concentration, and the growth factor VEGF diffusion reaction equation of the retinal unit grid, and save the information of the corresponding unit grid;

[0168] The diffusion reaction equation of the growth factor VEGF in the embodiment of the present application is obtained by a mass balance including diffusion, accumulation, and consumption terms. The situation in which the endothelial cells consume the growth factor VEGF (receptor-mediated binding) only occurs when n(x, y, t) = 2 and the endothelial cells consume the growth factor VEGF at a rate of λ, and is proportional to λ. When the cell grid attribute is ambient (n = 0), the growth factor VEGF is generated according to the oxygen concentration; specifically, as follows:

[0169]

[0170] in:

[0171] c υ is the current growth factor VEGF concentration of the retinal unit grid; the initial value is the steady-state value under normal conditions;

[0172] D is the diffusion coefficient matrix of growth factor VEGF, specifically:

[0173]

[0174] d is a fixed value;

[0175] λ * is the half-life degradation rate, which is a constant value;

[0176] λ is the consumption rate of growth factor VEGF by endothelial cells, which is a constant value;

[0177] η is the production efficiency of growth factor VEGF under the control of oxygen concentration, which is a constant value;

[0178] S(C o ) is the control equation for generating the growth factor concentration under the influence of the current oxygen concentration of the retinal unit grid; if the oxygen concentration of the unit grid is C lim The above does not produce growth factor VEGF due to lack of oxygen. If the oxygen concentration of the unit grid is less than C lim VEGF will be produced due to lack of oxygen, so the equation is as follows:

[0179]

[0180] Step S400: predict the movement direction of the blood vessel tip in the two-dimensional discrete coordinate system within the current time period t based on the fibronectin concentration of the retinal unit grid after t time, the growth factor concentration after t time, and the control equation of the movement direction of the blood vessel tip, and save it under the information of the corresponding unit grid.

[0181] The dynamics of the vascular tip are influenced by the structure and chemical properties of the extracellular matrix. The direction of movement of the vascular tip u is related to the growth factor concentration c after time t. υ and fibronectin concentration c after time t f The gradient of K is proportional to the properties of the extracellular matrix defined by K.

[0182] According to the gradient of growth factor VEGF concentration and fibronectin concentration, the control equation of the direction of blood vessel tip movement is:

[0183]

[0184] in:

[0185] is the chemotactic parameter of growth factor VEGF;

[0186] is the fibrin chemotaxis parameter;

[0187] K is the extracellular matrix anisotropy tensor matrix, which is an intrinsic property of the matrix space of each extracellular unit grid (x, y)∈Ω and is independent of time. It is non-uniform and anisotropic in space. This property means that the extension direction of the tip cell (blood vessel tip) is not necessarily consistent with the direction of the growth factor concentration gradient or the fibronectin concentration gradient; the governing equation is:

[0188]

[0189] k a is the anisotropy parameter, representing the resistance to movement along a certain direction;

[0190] k cond Heterogeneity parameters;

[0191] υ x and υ y is the component of the unit vector υ that describes the orientation of the extracellular matrix at the unit grid (x, y). At each unit grid (x, y), is a regular equation, randomly generating υ x and υ y , and is fixed at all times t. Similarly, for each position (x, y) ∈ Ω, k cond and k a Also randomly generated.

[0192] Step S500: Predict the first velocity of each blood vessel tip in the two-dimensional discrete coordinate system when moving in the first direction of motion based on the governing equations for the fibronectin concentration at the blood vessel tip t time later, the first direction of motion, and the velocity of the blood vessel tip. Predict the second velocity of each blood vessel tip in the two-dimensional discrete coordinate system when moving in the second direction of motion based on the governing equations for the fibronectin concentration at the blood vessel tip t time later, the second direction of motion, and the velocity of the blood vessel tip, and save the information in the corresponding unit grid; wherein the first direction of motion is the direction of motion of the blood vessel tip during the current time period t; and the second direction of motion is the direction of motion of the blood vessel tip during the previous time period t. In this embodiment of the present application, the governing equations for the velocity of the blood vessel tip are as follows:

[0193]

[0194] in:

[0195] P is the velocity of the blood vessel tip;

[0196] kp is the rate of vascular tip division, which is related to the concentration of growth factor VEGF c υ Related.

[0197]

[0198] l e is a constant, the length of the endothelial cell;

[0199] t c The cell division cycle is determined by the concentration of growth factor VEGF c υ Determine that the control equation is as follows:

[0200]

[0201] τ is a constant, the cell proliferation time coefficient;

[0202] is the average growth factor (VEGF) concentration; the sum of the growth factor concentrations of all unit grids in the two-dimensional discrete coordinate system at the current time is divided by the total number of unit grids, and is calculated once for each time period t;

[0203] c υ is the current concentration of growth factor VEGF;

[0204] c υlim is the threshold concentration of growth factor. When the VEGF concentration is higher than this value, the equation begins to work. It represents the critical time t c Depends on the current VEGF concentration c υ and average concentration When the VEGF concentration c υ Greater than the threshold c υlim When the critical time t c It will be calculated according to this equation.

[0205] Step S600: Calculate a first predicted extension distance of each blood vessel tip in a two-dimensional discrete coordinate system in a first direction of motion based on the discrete coordinates (x0, y0) of the corresponding unit grid of the blood vessel tip and the first motion velocity, and calculate a second predicted extension distance of each blood vessel tip in a second direction of motion based on the discrete coordinates (x0, y0) of the corresponding unit grid of the blood vessel tip and the second motion velocity, and save the calculations in the corresponding unit grid information; wherein the predicted extension distance is the total number of unit grids that the blood vessel tip passes through when it starts moving from the discrete coordinates (x0, y0) of the blood vessel tip in the corresponding direction of motion at the corresponding motion velocity and reaches the end coordinates (x1, y1) after t time, under the premise that no blood vessel branching or fusion occurs during extension. In this embodiment of the present application, step S600 specifically includes:

[0206] According to the movement speed of the blood vessel tip (which can be the first movement speed or the second movement speed), the vector (x p ,y p ) and the discrete coordinates (x0, y0) of the cell grid corresponding to the blood vessel tip to calculate the end point coordinates (x1, y1) of the cell grid this extension; where x1 = X0 + X p , y1=y0+y p ;

[0207] The Bresenham algorithm is used to calculate the discrete coordinates of the unit grids passing through the straight line connecting the discrete coordinates (x0, y0) of the unit grid corresponding to the blood vessel tip to the end point coordinates (x1, y1) of the current extension of the blood vessel tip, and these passed discrete coordinates of the unit grids are used as the predicted extension distance of the blood vessel tip (x0, y0) in the corresponding movement direction (the first predicted extension distance when the movement speed is the first movement speed, and the second predicted extension distance when the movement speed is the second movement speed).

[0208] Step S700: Perform an extension simulation for each blood vessel tip in the two-dimensional discrete coordinate system within a current time period t based on the first movement direction of the blood vessel tip, the existence time of the blood vessel tip, the second movement direction of the blood vessel tip, the first movement speed, the second movement speed, the first predicted extension distance, and the second predicted extension distance. The unit cell attributes and the existence time of each unit cell in the two-dimensional discrete coordinate system are updated as the extension simulation of each blood vessel tip is performed. The simulation time is updated when the extension simulation of all blood vessel tips in the two-dimensional discrete coordinate system is completed. The extension simulation includes blood vessel extension, blood vessel branching that occurs during blood vessel extension, and blood vessel fusion that occurs during blood vessel extension.

[0209] Blood vessel branching occurs when a single blood vessel tip splits into two new ones. This is an important mechanism for the formation of new blood vessels during angiogenesis. When branching occurs, a single blood vessel tip becomes two new ones. In this example, to simulate this process, the following definitions are defined:

[0210] (1) The conditions for vascular branching are as follows:

[0211] 1. The existence time of the blood vessel tip is greater than a preset blood vessel tip existence time threshold τ a ;

[0212] 2. Unit velocity vector of the blood vessel tip during this time period t The angle change between the angle between the positive (or negative) direction of the x-axis of the coordinate system and the unit velocity vector of the blood vessel tip (i.e., the original motion direction of the blood vessel tip) and the positive (or negative) direction of the x-axis of the coordinate system in the last time period t is greater than the preset angle threshold τ u .

[0213] (2) When both conditions for vascular branching are met simultaneously, vascular branching at the vascular tip is simulated as follows:

[0214] A new vascular tip is generated from an existing vascular tip, that is, the same unit grid with the vascular tip attribute has two vascular tip extension attributes, that is, the movement direction of the vascular tip corresponding to the unit grid in the current time period t is set to the first movement direction and the second movement direction, and starting from the vascular tip where the vascular branch occurs, it is extended in the original movement direction (second movement direction) according to the second predicted extension distance; and it is extended in the newly predicted movement direction (first movement direction) according to the first predicted extension distance.

[0215] Vascular fusion occurs when two vascular tips, or one vascular tip, meet and fuse with another. In tip-to-tip fusion, both or only one of the vascular tips may disappear, depending on the specific circumstances. In tip-to-tip fusion, the vascular tip disappears as part of a ring, no longer extending, meaning the tip no longer exists. To simulate this process, the following definitions are defined in this embodiment:

[0216] (1) The conditions for vascular fusion are as follows:

[0217] Among the unit grids included in the predicted extension distance of the blood vessel tip (including the first predicted extension distance and / or the second predicted extension distance), there are unit grids whose unit grid attributes are blood vessels and / or blood vessel tips.

[0218] (2) When the vascular fusion conditions are met, vascular fusion occurs at the simulated vascular tip as follows:

[0219] In the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point on the corresponding predicted extension distance (the order of the unit grids and unit grid attributes on the predicted extension distance and the actual extension distance mentioned in this embodiment is in this order), the blood vessel tip is sequentially fused with the blood vessel and the blood vessel tip on the predicted extension distance, and extended according to the actual extension distance of the blood vessel tip after the blood vessel fusion, until the extension simulation of the blood vessel tip in this time period t ends. The actual extension distance is all unit grids (including the unit grid (x0, y0) to which the blood vessel tip belongs) passes from the discrete coordinates (x0, y0) of the blood vessel tip to the end of the extension simulation after blood vessel branching and / or blood vessel fusion occurs during the process of the blood vessel tip starting from the discrete coordinates (x0, y0) of the blood vessel tip in the corresponding movement direction and at the corresponding movement speed according to the corresponding predicted extension distance. In this embodiment of the present application, the blood vessel tip is fused with the blood vessel and the blood vessel tip on the predicted extension distance, and extended according to the actual extension distance of the blood vessel tip after the blood vessel fusion includes:

[0220] 1. If the vessel fused with the vessel tip is within the predicted extension distance (the first predicted extension distance and / or the second predicted extension distance), it means that the vessel tip needs to disappear at the unit grid where the vessel is located, and the extension of the vessel tip will also stop at the unit grid where the vessel is located. That is, the actual extension distance of the vessel tip is set to the unit grid to which the vessel tip belongs and all unit grids before the unit grid to which the vessel tip fused with the vessel tip within the predicted extension distance (including the unit grid to which the vessel belongs). The vessel tip existence time of the vessel tip is reset to zero and the timing is stopped. The unit grid attributes of all unit grids within the actual extension distance are changed to vessels, and the extension simulation of the vessel tip within this time period t is terminated.

[0221] 2. If the vessel tip that undergoes vascular fusion is the vessel tip at the predicted extension distance (the first predicted extension distance and / or the second predicted extension distance), vascular fusion is performed based on the movement direction of the vessel tip and the movement direction of the vessel tip at the predicted extension distance, specifically:

[0222] (1) If the angle between the movement direction of the vascular tip (which can be the first movement direction or the second movement direction) and the movement direction of the vascular tip at the corresponding predicted extension distance (which can be the first predicted extension distance or the second predicted extension distance) is greater than or equal to 90 degrees, it is considered that both tips become inactive. At this time, the vascular tip where the vascular fusion occurs and the vascular tip at its predicted extension distance will disappear, that is, the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids before the unit grid to which the vascular tip belongs where the vascular fusion occurs at the predicted extension distance. The vascular tip existence time of the vascular tip and the vascular tip at the predicted extension distance is cleared and the timing is stopped. The unit grid attributes of all unit grids at the actual extension distance are changed to blood vessels, and the extension simulation of the vascular tip in this time period t is ended.

[0223] (2) If the angle between the movement direction of the vascular tip and the movement direction of the vascular tip at the corresponding predicted extension distance is less than 90 degrees, it is considered that only one of the two vascular tips that have fused has become inactive. One vascular tip is randomly selected to continue extending at the corresponding predicted extension distance, and the other vascular tip stops extending at the unit grid where the two vascular tips have fused. Since the two tips are in the same position in the model (i.e., the unit grid where the two vascular tips have fused), the choice of which tip becomes inactive will not affect the determinism of the model.

[0224] If the vascular tip is randomly selected for extension, the vascular tip on the predicted extension distance stops extending; when there is no more blood vessel or vascular tip after the vascular tip on the predicted extension distance that merges with the vascular tip on the corresponding predicted extension distance, the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids included in the predicted extension distance of the vascular tip, and the vascular tip existence time of the vascular tip and the vascular tip on the predicted extension distance that merges with it is cleared and the timing is stopped, the unit grid attribute of the last unit grid on the actual extension distance is changed to the vascular tip and the timer is started to record the vascular tip existence time from the initial value and the other unit grids on the actual extension distance are reset. The cell mesh attributes of the remaining cell meshes are changed to blood vessels; the extension simulation of the blood vessel tip within the current time period t is terminated; when there are still blood vessels and / or blood vessel tips after the blood vessel tip at the predicted extension distance that is fused with the blood vessel tip at the corresponding predicted extension distance, the blood vessel tip continues to fuse with the blood vessel and blood vessel tip at the predicted extension distance, and extends according to the actual extension distance of the blood vessel tip after blood vessel fusion, until the extension simulation of the blood vessel tip within the current time period t is terminated; wherein, the initial value of the timer of the last cell mesh at the actual extension distance is the existence time of the blood vessel tip before the blood vessel tip is cleared if no blood vessel branching occurs, and the initial value is zero if a blood vessel branching occurs;

[0225] If the vascular tip is randomly selected to stop extending, the vascular tip is extended at the predicted extension distance where the vascular tip fuses; the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids before the unit grid to which the vascular tip belongs at the predicted extension distance where the vascular tip fuses (excluding the unit grid to which the vascular tip belongs at the predicted extension distance), and the vascular tip existence time of the vascular tip is cleared and the timing is stopped, and the attributes of all unit grids at the actual extension distance are changed to blood vessels, ending the extension simulation of the vascular tip within this t time period; regardless of whether the vascular tip at the predicted extension distance where the vascular tip fuses has completed the extension simulation within this t time period, the vascular tip extension simulation continues, and it has no effect on it.

[0226] In the embodiment of the present application, in order to simulate the process of blood vessel extension, the following definition is given:

[0227] Each blood vessel tip in the two-dimensional discrete coordinate system is sequentially extended within the current time period t by simulation. The extension simulation of the blood vessel tip within the current time period t includes:

[0228] Determine whether vascular branching occurs at the vascular tip;

[0229] If no branching occurs at the vascular tip, it is determined whether vascular fusion occurs at the vascular tip within the first predicted extension distance. If no vascular fusion occurs, the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids within the first predicted extension distance. The current vascular tip existence time of the vascular tip is used as the initial value of the timer for the last unit grid within the actual extension distance. The vascular tip existence time of the vascular tip is reset and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to vascular tip, and the timer is started to record the vascular tip existence time from the initial value. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessels, and the extension simulation of the vascular tip within the current time period t is terminated. If vascular fusion occurs, the vascular tip is sequentially fused with the blood vessels and the vascular tip within the first predicted extension distance in the order from the discrete coordinates of the vascular tip to the end point coordinates within the first predicted extension distance. The vascular tip is then extended according to the actual extension distance of the vascular tip after vascular fusion, until the extension simulation of the vascular tip within the current time period t is terminated.

[0230] If a blood vessel branch occurs at the tip of the blood vessel, the tip of the blood vessel is extended in the second movement direction according to the second predicted extension distance and the tip of the blood vessel is extended in the first movement direction according to the first predicted extension distance; wherein:

[0231] The extending of the blood vessel tip in the first movement direction according to the first predicted extension distance specifically includes:

[0232] determining whether vascular fusion occurs at the tip of the blood vessel at the first predicted extension distance;

[0233] If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids in the first predicted extension distance, the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped, the unit grid attribute of the last unit grid in the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero, and the unit grid attributes of the remaining unit grids in the actual extension distance are changed to blood vessel, and the extension simulation of the blood vessel tip in this time period t is ended;

[0234] If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips at the first predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point at the first predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed;

[0235] The extending of the second movement direction of the blood vessel tip according to the second predicted extension distance specifically includes:

[0236] determining whether vascular fusion occurs at the tip of the blood vessel at the second predicted extension distance;

[0237] If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids within the second predicted extension distance. The blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessels, and the extension simulation of the blood vessel tip within the current time period t is terminated.

[0238] If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips on the second predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point on the second predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed.

[0239] like Figure 2 As shown, in this embodiment, step S700 includes the following steps:

[0240] Step S701: Determine whether a blood vessel branch occurs during the extension simulation of the blood vessel tip within this time period t based on the first movement direction of the blood vessel tip, the existence time of the blood vessel tip, the second movement direction of the blood vessel tip, and the blood vessel branching condition. If a blood vessel branch occurs, execute step S711; otherwise, execute step S721.

[0241] Step S711, the blood vessel tip is extended according to the first predicted extension distance and the second predicted extension distance respectively, and then steps S712 and S713 are executed;

[0242] Step S712: Determine whether vascular fusion occurs when the vascular tip extends according to the first predicted extension distance based on the unit grid attributes of each unit grid and the vascular fusion condition; if vascular fusion occurs, execute step S714; otherwise, execute step S715;

[0243] Step S713: determining whether vascular fusion occurs when extending according to the second predicted extension distance based on the unit grid attributes of each unit grid and the vascular fusion condition; if vascular fusion occurs, executing step S716; otherwise, executing step S717;

[0244] Step S714: The blood vessel tip is sequentially fused with the blood vessels and blood vessel tips within the first predicted extension distance in the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point within the first predicted extension distance. The blood vessel tip is then extended according to the actual extension distance of the blood vessel tip after fusion, until the extension simulation of the blood vessel tip within the current time period t is completed. Then, step S702 is executed.

[0245] Step S715: The actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids within the first predicted extension distance. The blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessel, and the extension simulation of the blood vessel tip within the current time period t is terminated. Then, step S702 is executed.

[0246] Step S716: The blood vessel tip is sequentially fused with the blood vessels and blood vessel tips at the second predicted extension distance in the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point at the second predicted extension distance. The blood vessel tip is extended according to the actual extension distance of the blood vessel tip after the fusion, until the extension simulation of the blood vessel tip in the current time period t is completed; then step S702 is executed.

[0247] Step S717: The actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids within the second predicted extension distance. The blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessel, and the extension simulation of the blood vessel tip within the current time period t is terminated. Then, step S702 is executed.

[0248] Step S721: Determine whether vascular fusion occurs at the tip of the blood vessel within the first predicted extension distance. If vascular fusion does not occur, execute step S722; otherwise, execute step S723.

[0249] Step S722: The blood vessel tip is sequentially fused with the blood vessels and blood vessel tips within the first predicted extension distance in the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point within the first predicted extension distance. The blood vessel tip is then extended according to the actual extension distance of the blood vessel tip after fusion, until the extension simulation of the blood vessel tip within the current time period t is completed. Then, step S702 is executed.

[0250] Step S723: The actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids within the first predicted extension distance. The current blood vessel tip existence time of the blood vessel tip is used as the initial value of the timer for the last unit grid within the actual extension distance. The blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from the initial value. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessel. The extension simulation of the blood vessel tip within the current time period t is terminated, and then step S702 is executed.

[0251] Step S702: Determine whether all vascular tips in the two-dimensional discrete coordinate system have completed the extension simulation within the current time period t. If so, execute step S703; otherwise, extend the next vascular tip in the two-dimensional discrete coordinate system within the current time period t and execute step S701.

[0252] Step S703 , the simulation time T is added to the preset time t, and the unit grid properties of each unit grid in the two-dimensional discrete coordinate system after the time t are saved, and then step S704 is executed.

[0253] Step S704: Determine whether the simulation termination condition has been met. If the simulation termination condition has been met, output the simulation results and terminate. Otherwise, enter the next extended simulation for a time period t, and repeat the above process starting from step S300; wherein the simulation results include but are not limited to the simulation time T, the initial value of the patient's blood glucose concentration, the initial value of the oxygen concentration, the unit grid attributes of each unit grid in the two-dimensional discrete coordinate system after the last extension within the time period t, and the unit grid attributes of each unit grid in the two-dimensional discrete coordinate system after each extension within the time period t. The simulation termination condition can be a preset simulation time threshold, or a simulation time T input by the user according to needs, or the percentage of unit grids in the two-dimensional discrete coordinate system that are blood vessels to the total number of unit grids in the two-dimensional discrete coordinate system, and the embodiments of the present application do not limit this.

[0254] The embodiment of the present application provides a method for simulating the diabetic retinopathy process, which fully considers the impact of diabetic patients' blood sugar levels on diabetic retinopathy and incorporates oxygen concentration and blood sugar concentration into the simulated vascular extension model. It not only simulates the occurrence and development process of diabetic retinopathy more comprehensively and understands how these factors interact to cause retinal vascular abnormalities, which helps to reveal the root cause of the disease, but also achieves accurate simulation prediction of diabetic retinopathy, providing a reliable basis for optimizing and adjusting diabetic retinopathy treatment plans and developing and verifying new drugs. In addition, the embodiment of the present application also considers the blood sugar control level and oxygen supply status of individual patients, establishes personalized simulation models for different patients, and realizes personalized disease prediction and treatment plan optimization.

[0255] Obviously, those skilled in the art may make various changes and modifications to the present invention without departing from the spirit and scope of the present invention. Thus, if such changes and modifications fall within the scope of the claims and their equivalents, the present invention is intended to include such changes and modifications.

Claims

1. A method for simulating the process of diabetic retinopathy, characterized in that: The following steps are involved: Step S100: constructing a two-dimensional discrete coordinate system and defining a blood vessel index function n(x, y, t) that describes and records the cell grid attributes of each cell grid in the two-dimensional discrete coordinate system at each preset time t; wherein the cell grid attributes include the blood vessel, the blood vessel tip, and the environment; Step S200: Initializing the blood glucose concentration and oxygen concentration of each unit grid in the two-dimensional discrete coordinate system according to the patient's examination results, and initializing the unit grid attributes, fibronectin concentration, and growth factor concentration of each unit grid in the two-dimensional discrete coordinate system according to a preset initial state; Step S300: updating the blood glucose concentration, oxygen concentration, fibronectin concentration, and growth factor concentration of each unit grid at a preset time interval t based on the unit grid's current blood glucose concentration, the governing equation for the blood glucose concentration, the current oxygen concentration, the governing equation for the oxygen concentration, the current unit grid attributes, the current fibronectin concentration, the governing equation for the fibronectin concentration, the current growth factor concentration, and the growth factor diffusion reaction equation; Step S400, predicting the movement direction of each blood vessel tip in the two-dimensional discrete coordinate system within the current time period t based on the fibronectin concentration of the unit grid after time t, the growth factor concentration after time t, and the blood vessel tip movement direction control equation; Step S500: Predicting a first motion velocity of each blood vessel tip in a two-dimensional discrete coordinate system when moving in a first motion direction based on a governing equation for the fibronectin concentration at the blood vessel tip after time t, the first motion direction, and the motion velocity of the blood vessel tip; and predicting a second motion velocity of each blood vessel tip in a two-dimensional discrete coordinate system when moving in a second motion direction based on a governing equation for the fibronectin concentration at the blood vessel tip after time t, the second motion direction, and the motion velocity of the blood vessel tip; wherein the first motion direction is the motion direction of the blood vessel tip during the current time period t; and the second motion direction is the motion direction of the blood vessel tip during the previous time period t. Step S600, calculating a first predicted extension distance of each blood vessel tip in a two-dimensional discrete coordinate system in a first motion direction based on the discrete coordinates (x0, y0) of the blood vessel tip and a first motion speed, and calculating a second predicted extension distance of each blood vessel tip in a second motion direction based on the discrete coordinates (x0, y0) of the blood vessel tip and a second motion speed; wherein the first predicted extension distance or the second predicted extension distance is the distance that the blood vessel tip moves in the first motion direction or the second motion direction at the first motion speed or the second motion speed starting from the discrete coordinates (x0, y0) of the unit grid to which the blood vessel tip belongs, and reaching the end coordinate after a preset time t. All cell grids passed; Step S700: Perform an extension simulation on each blood vessel tip in the two-dimensional discrete coordinate system within the current time period t based on the first movement direction of the blood vessel tip, the existence time of the blood vessel tip, the second movement direction of the blood vessel tip, the first movement speed, the second movement speed, the first predicted extension distance, and the second predicted extension distance. The unit grid attributes and the existence time of each unit grid in the two-dimensional discrete coordinate system are updated as the extension simulation of each blood vessel tip is performed. The simulation time is updated when the extension simulation of all blood vessel tips in the two-dimensional discrete coordinate system is completed. The extension simulation includes blood vessel extension, blood vessel branching that occurs during blood vessel extension, and blood vessel fusion that occurs during blood vessel extension.

2. The method according to claim 1, characterized in that The step S300 includes the following steps: Predict the blood glucose concentration of each unit grid after time t based on the current blood glucose concentration of the unit grid and the control equation of the blood glucose concentration; and save the information under the corresponding unit grid; Predict the oxygen concentration of each unit grid after time t based on the current oxygen concentration of the retinal unit grid, the blood glucose concentration after time t, and the control equation of the oxygen concentration; and save it under the information of the corresponding unit grid; Predicting the fibronectin concentration of each retinal unit grid after time t based on the current unit grid properties of the retinal unit grid, the current fibronectin concentration and the control equation of the fibronectin concentration, and saving it under the information of the corresponding unit grid; According to the current oxygen concentration, growth factor concentration and growth factor diffusion reaction equation of the unit grid, the growth factor concentration of each unit grid after time t is predicted, and the information of the corresponding unit grid is saved.

3. The method according to claim 1, characterized in that The governing equation for the oxygen concentration is as follows: in: is the current oxygen concentration of the retinal cell grid; is the diffusion coefficient of oxygen in retinal blood vessels, which is a fixed constant; The oxygen consumption rate equation for the retinal cell grid is: A is the normal oxygen consumption rate of the unit grid, which is a fixed constant; B is the constant of oxygen consumption of the unit grid under the influence of blood glucose; The lowest oxygen concentration for normal physiological activities of the unit grid is a fixed constant. Unit grids below this constant value are in an oxygen-deficient state. is the amount of oxygen transported from the surrounding blood vessels in the unit grid. Vascular transport is the source of oxygen for the retinal unit grid. According to the current blood vessel density, The calculation equation is as follows: is the oxygen transport rate constant in blood vessels; C is the transport efficiency constant under the influence of blood glucose; The retinal cell grid At the current time The blood vessel density at the time, specifically the current unit grid and the number of cell grids containing blood vessels in the nine cell grids adjacent to the current cell grid.

4. The method according to claim 1, wherein The governing equation for the fibronectin concentration is: in: is the current fibronectin concentration of the retinal unit grid; the initial value has been set in step S200; is the production rate of fibronectin in endothelial cells, which is a fixed constant; is the uptake rate of fibronectin by endothelial cells, which is a fixed constant.

5. The method according to claim 1, wherein The growth factor diffusion reaction equation is as follows: in: is the current growth factor VEGF concentration of the retinal unit grid; the initial value is the steady-state value under normal conditions; is the diffusion coefficient matrix of growth factor VEGF, specifically: d is a fixed value; is the half-life degradation rate, which is a constant value; is the consumption rate of growth factor VEGF by endothelial cells, which is a constant value; is the production efficiency of growth factor VEGF under the control of oxygen concentration, which is a constant value; To generate the governing equation for the growth factor concentration under the influence of the current oxygen concentration of the retinal cell grid; if the oxygen concentration of the cell grid is The above does not produce growth factor VEGF due to lack of oxygen. If the oxygen concentration of the cell grid is less than VEGF will be produced due to lack of oxygen, so the equation is as follows: in, It is the minimum oxygen concentration for normal physiological activities of the unit grid, which is a fixed constant. The unit grid below this constant value is in an oxygen-deficient state.

6. The method according to claim 1, characterized in that The control equation for the direction of movement of the blood vessel tip is as follows: in: is the current growth factor VEGF concentration of the retinal unit grid; the initial value is the steady-state value under normal conditions; is the current fibronectin concentration of the retinal unit grid; the initial value has been set in step S200; is the chemotactic parameter of growth factor VEGF; is the fibrin chemotaxis parameter; is the extracellular matrix anisotropy tensor matrix, which is an intrinsic property of the unit grid (x, y) outside each cell and is independent of time, where (x, y)∈Ω; it is non-uniform and anisotropic in space. This property means that the extension direction of the tip cell is not necessarily consistent with the direction of the growth factor concentration gradient or the fibronectin concentration gradient; the governing equation is: is the anisotropy parameter, representing the resistance to movement along a certain direction; Heterogeneity parameters; and is the description unit grid The unit vector of the extracellular matrix orientation at Component, in each unit grid ,by is a regular equation, randomly generated and , and at all times Fixed; for each unit grid ∈ Ω, and Also randomly generated.

7. The method according to claim 1, characterized in that The control equation of the blood vessel tip movement speed is as follows: in: u is the direction of movement of the blood vessel tip; P is the velocity of the blood vessel tip; is the rate of vascular tip division and is related to the concentration of growth factor VEGF Related; is a constant, the length of the endothelial cell; For the cell division cycle, the concentration of growth factor VEGF Determine that the control equation is as follows: is a constant, the cell proliferation time coefficient; is the average growth factor (VEGF) concentration; the sum of the growth factor concentrations of all unit grids in the two-dimensional discrete coordinate system at the current time is divided by the total number of unit grids, and is calculated once for each time period t; is the current concentration of growth factor VEGF; is the threshold concentration of growth factors.

8. The method according to claim 1, characterized in that The step S600 includes the following steps: The vectors in the x and y directions according to the velocity of the blood vessel tip and the discrete coordinates of the vessel tip Calculate the coordinates of the end point of this extension of the cell grid ; where x1=x0+x p , y1=y0+y p ; The movement speed is the first movement speed or the second movement speed; Use Bresenham algorithm to calculate the discrete coordinates of the blood vessel tip The coordinates of the end point of this extension to the blood vessel tip The discrete coordinates of the unit grids passed by this straight line are used as the vascular tip. The predicted extension distance within this time period t; the predicted extension distance is the first predicted extension distance or the second predicted extension distance.

9. The method according to claim 8, characterized in that The conditions for the vascular tip to branch when extending include: The vessel existence time of the vessel tip is greater than a preset vessel tip existence time threshold ; The unit velocity vector of the blood vessel tip during this time period t The angle change between the angle between the positive direction of the x-axis of the coordinate system and the angle between the unit velocity vector of the blood vessel tip and the positive direction of the x-axis of the coordinate system in the last time period t is greater than the preset angle threshold ; When the two conditions for vascular branching are met at the same time, the process of vascular branching at the vascular tip is simulated as follows: Starting from the tip of the blood vessel where the blood vessel branch occurs, the blood vessel is extended in the second movement direction according to the second predicted extension distance and in the first movement direction according to the first predicted extension distance; The conditions for vascular fusion to occur when the vascular tip extends include: Among the cell grids included in the blood vessel tip prediction extension distance, there are cell grids whose cell grid attributes are blood vessels and / or blood vessel tips; When the vascular fusion conditions are met, the vascular fusion process at the vascular tip is simulated including: According to the order from the discrete coordinates of the blood vessel tip to the coordinates of the end point on the corresponding predicted extension distance, the blood vessel tip is sequentially fused with the blood vessels and blood vessel tips on the predicted extension distance, and the blood vessel tip is extended according to the actual extension distance after the blood vessel fusion, until the extension simulation of the blood vessel tip in this time period t is completed; The blood vessel tip extending the blood vessel comprises: Each blood vessel tip in the two-dimensional discrete coordinate system is sequentially extended within the current time period t by simulation, wherein the blood vessel tip extension simulation within the current time period t is: Determine whether vascular branching occurs at the vascular tip; If no branching occurs at the vascular tip, it is determined whether vascular fusion occurs at the vascular tip within the first predicted extension distance. If no vascular fusion occurs, the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids within the first predicted extension distance. The current vascular tip existence time of the vascular tip is used as the initial value of the timer for the last unit grid within the actual extension distance. The vascular tip existence time of the vascular tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid within the actual extension distance is changed to vascular tip, and the timer is started to record the vascular tip existence time from the initial value. The unit grid attributes of the remaining unit grids within the actual extension distance are changed to blood vessels, and the extension simulation of the vascular tip within the current time period t is terminated. If vascular fusion occurs, the vascular tip is sequentially fused with the blood vessels and the vascular tip within the first predicted extension distance in the order from the discrete coordinates of the vascular tip to the end point coordinates within the first predicted extension distance. The vascular tip is extended according to the actual extension distance of the vascular tip after vascular fusion until the extension simulation of the vascular tip within the current time period t is terminated. If a blood vessel branch occurs at the tip of the blood vessel, the tip of the blood vessel is extended in the second movement direction according to the second predicted extension distance and the tip of the blood vessel is extended in the first movement direction according to the first predicted extension distance; wherein: The extending of the blood vessel tip in the first movement direction according to the first predicted extension distance specifically includes: determining whether vascular fusion occurs at the tip of the blood vessel at the first predicted extension distance; If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids in the first predicted extension distance, the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped, the unit grid attribute of the last unit grid in the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero, and the unit grid attributes of the remaining unit grids in the actual extension distance are changed to blood vessel, and the extension simulation of the blood vessel tip in this time period t is ended; If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips at the first predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point at the first predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed; The extending of the second movement direction of the blood vessel tip according to the second predicted extension distance specifically includes: determining whether vascular fusion occurs at the tip of the blood vessel at the second predicted extension distance; If blood vessel fusion does not occur, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids on the second predicted extension distance, and the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped. The unit grid attribute of the last unit grid on the actual extension distance is changed to blood vessel tip, and the timer is started to record the blood vessel tip existence time from zero. The unit grid attributes of the remaining unit grids on the actual extension distance are changed to blood vessels, and the extension simulation of the blood vessel tip in this time period t is ended. If vessel fusion occurs, the vessel tip is sequentially fused with the vessels and vessel tips on the second predicted extension distance in the order from the discrete coordinates of the vessel tip to the coordinates of the end point on the second predicted extension distance, and the vessel tip is extended according to the actual extension distance of the vessel tip after vessel fusion until the extension simulation of the vessel tip in this time period t is completed.

10. The method according to claim 9, characterized in that The vascular tip is fused with the blood vessel and the vascular tip at the predicted extension distance, and extended according to the actual extension distance of the vascular tip after the vascular fusion, including: If the vessel fused with the blood vessel tip is within the predicted extension distance, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids before the unit grid to which the blood vessel tip belongs within the predicted extension distance, the blood vessel tip existence time of the blood vessel tip is reset to zero and the timer is stopped, and the unit grid attributes of all unit grids within the actual extension distance are changed to blood vessels, thereby ending the extension simulation of the blood vessel tip within the current time period t. If the blood vessel tip that undergoes vascular fusion is the blood vessel tip at the predicted extension distance, the blood vessel fusion is performed according to the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the predicted extension distance, specifically: If the angle between the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the corresponding predicted extension distance is greater than or equal to 90 degrees, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids before the unit grid to which the blood vessel fusion occurs at the predicted extension distance. The blood vessel tip existence time of the blood vessel tip and the blood vessel tip at the predicted extension distance is reset to zero and the timing is stopped. The unit grid attributes of all unit grids at the actual extension distance are changed to blood vessels, and the extension simulation of the blood vessel tip in this time period t is terminated. If the angle between the movement direction of the blood vessel tip and the movement direction of the blood vessel tip at the corresponding predicted extension distance is less than 90 degrees, if the blood vessel tip is randomly selected for extension, the blood vessel tip at the predicted extension distance stops extending; when there is no blood vessel or blood vessel tip after the blood vessel tip at the predicted extension distance that merges with the blood vessel tip at the corresponding predicted extension distance, the actual extension distance of the blood vessel tip is set to the unit grid to which the blood vessel tip belongs and all unit grids included in the predicted extension distance of the blood vessel tip, and the blood vessel tip existence time of the blood vessel tip and the blood vessel tip at the predicted extension distance that merges with it is cleared and the timing is stopped, the unit grid attribute of the last unit grid at the actual extension distance is changed to the blood vessel tip and the timer is started to count from the initial value Recording the existence time of the blood vessel tip and changing the unit grid attributes of the remaining unit grids in the actual extension distance to blood vessels, ending the extension simulation of the blood vessel tip in this time period t; when there are still blood vessels and / or blood vessel tips after the blood vessel tip in the predicted extension distance that merges with the blood vessel tip in the corresponding predicted extension distance, the blood vessel tip continues to fuse with the blood vessel and blood vessel tip in the predicted extension distance, and extends according to the actual extension distance of the blood vessel tip after blood vessel fusion, until the extension simulation of the blood vessel tip in this time period t is ended; wherein, the initial value of the timer of the last unit grid in the actual extension distance is the existence time of the blood vessel tip before the blood vessel tip is cleared if no blood vessel branch occurs, and the initial value is zero if a blood vessel branch occurs; If the vascular tip is randomly selected to stop extending, the vascular tip is extended at the predicted extension distance where the vascular tip merges with the vascular tip; the actual extension distance of the vascular tip is set to the unit grid to which the vascular tip belongs and all unit grids before the unit grid to which the vascular tip belongs at the predicted extension distance where the vascular tip merges with the vascular tip at the corresponding predicted extension distance, and the vascular tip existence time of the vascular tip is cleared and the timing is stopped, and the attributes of all unit grids at the actual extension distance are changed to blood vessels, ending the extension simulation of the vascular tip within this time period t.

Images