Evaluation of myocardial perfusion using non-invasive measurements and seepage mechanics

Abstract

A method of describing fluid dynamics in a tissue or organ using a numerical model of blood flow is provided. In one aspect, the invention relates to a method for assessing myocardial perfusion in a human patient using non-invasive measurements and percolation mechanics including computational fluid dynamics. The method may be implemented as a computer-implemented invention. The method includes using patient-specific anatomical structure data obtained from the non-invasive measurements to generate an anatomical geometric model of the myocardium connected to the plurality of coronary vessels. The anatomical geometric model is a combination of fluid and porous domains and boundaries thereof. The invention eliminates the deficiencies of known in vivo methods, including all risks associated with invasive procedures. Furthermore, an actually infinite temporal-spatial resolution is achieved, the result of the invention is completely free of motion artifacts, and the exposure of the patient to radiation is limited.

Description

Field of the Invention

[0001] The present invention relates to myocardial perfusion imaging (MPI) and permeability mechanics. In particular, the present invention relates to the use of non-invasive measurements (e.g., computed angiography or nuclear magnetic resonance) to construct a model of a porous domain of a human heart to be used for computational fluid dynamics. The results of the present invention allow diagnosis of the human heart. For example, diagnosis of coronary heart disease.

[0002] background

[0003] One of the main challenges for the application of perfusion assessment methods in medical diagnosis is coronary heart disease (CHD). According to the National Center for Health Statistics, 5.6% of American adults had CHD in 2018. The prevalence increases with age, and nearly 24% of people aged 75 years and older have CHD. Men are more likely to develop the disease than women - 7.4% and 4.1% of all American adults, respectively (Summary tables of health statistics for adults in the United States by Villarroel MA, Blackwell DL, Jen A (2019): 2018 National Health Interview Survey, National Center for Health Statistics, accessed on December 15, 2021). In 2018, CHD caused 365,744 deaths in the United States, almost 13% of all deaths (see Centers for Disease Control and Prevention, National Center for Health Statistics, National Vital Statistics System: public use data file documentation: mortality multiple cause-of-death micro-data files, accessed December 15, 2021). In the European Union, 623,276 people died from ischemic heart disease in 2017, almost 12% of all deaths in 2018 (see Eurostat 2021 Causes of death - deaths by country of residence and occurrence [HLTH_CD_ARO], accessed December 15, 2021).

[0004] The heart is the center of the circulatory system. The human heart beats about 100,000 times a day, a total of about 3 billion times in a lifetime. The pumping of the heart and the heartbeat are caused by the alternating contraction and relaxation of the myocardium. This pumping ensures the continuity of the movement of blood in the blood vessels. Therefore, the heart is responsible for the realization of various functions in the human body: supply of oxygen and nutrients, removal of metabolites, protection against foreign bodies, toxic substances and prevention of blood loss. In addition, the heart is responsible for maintaining constant physical and chemical conditions in the human body, especially temperature and pH (according to Ganong's Review of Medical Physiology, 26th edition, McGraw-Hill Education, Barrett KE et al., 2019). In fact, modern medical diagnostics provides a wide range of methods and parameters for the functional assessment of the heart. Nevertheless, just as more than a hundred years ago, a parameter that plays a key role in medical diagnosis, which describes the amount of blood pumped in a certain time interval (usually in a single cycle or one minute). As long as this parameter is maintained at a sufficient level, the prerequisites for other organs to maintain their functions are met. The recognition of this principle and attempts to evaluate this parameter can be found in very early physiological experimental work. As early as 1761, Albrecht von Haller studied pulmonary circulation by injecting colored fluids into the vena cava of fresh animal carcasses (see Haller A, Elementaphysiologiae corporis humani, Francisci Grasser & Sociorum, Lausanne, Switzerland, 1761).The turning point and foundation of modern blood flow diagnostic methods with indicator dilution technology (i.e., the technology commonly used today to assess perfusion) were the studies done by George Stewart (see Argueta EE, Paniagua D, 2019, Thermodilutioncardiac output, Cardiology in Review, 27(3), pp. 138-144), and then less than a century later, the studies done by Kenneth Zierler (see Meier P, Zierler KL, 1954, On the theory of the indicator-dilution method for measurement of blood flow and volume, The Journal of Applied Physiology, 6(12), pp. 731-744, Zierler K, 2000, Indicatordilution methods for measuring blood flow, volume, and other properties of biological systems: a brief history and memoir, Annals of Biomedical Engineering, 28(8), pp. 836-848).

[0005] Myocardial perfusion (ie, blood flow through a particular tissue mass) is a fundamental physiological parameter determined based on cardiology image data. Perfusion imaging, including myocardial perfusion imaging (MPI), involves injecting a contrast agent into the blood and recording its concentration profile over time in the cardiovascular system. A key issue in perfusion analysis is the calculation of so-called perfusion parameters or perfusion metrics such as myocardial blood flow MBF, myocardial blood volume MBV, mean transit time MTT and time to peak TTP (according to Taylor SH, Measurement of the cardiac output in man, Proc R Soc Med, Suppl 59 (Suppl 1), pp. 35-53, 1966, Axel L, Cerebral blood flow determination by rapid-sequence computed tomography: theoretical analysis, Radiology, 137 (3), pp. 679-86, 1980, Lee TY, Functional CT: Physiological models, Trends Biotechnol, 20, pp. S3-S10, Tomandl BF et al., Comprehensive imaging of ischemic stroke with multisection CT, Radiographics, 23 (3), pp. 65-92, Konstadt et al., 2003). AA et al., 2009, Theoretic basis and technical implementations of CT perfusion in acute ischemic stroke. Part 1: Theoretic basis, AJNR Am J Neuroradiol, 30(4), pp. 662-8). Let us assume that a certain control volume of interest (VOI) has a uniquely defined inlet and outlet. If a certain amount (m) of contrast agent is introduced at the initial time t0 = 0 and sampled at the outlet, then at time Δt i , after which the concentration of the contrast agent c(t i ),t i ∈Δt i , then the mass of contrast agent left in VOI is expressed by the following formula:

[0006] Δm(t i )=c(t i )ΔV(t i ) , (1)

[0007] Among them, ΔV(t i )=q(t i )·Δt i is the blood volume flowing through the outlet segment. After sampling at the outlet and determining the concentration over a sufficiently long period of time, the total mass of contrast agent that will reach the outlet is expressed by the following formula:

[0008]

[0009] If we take samples frequently enough, we get:

[0010]

[0011] Therefore, based on the mean value theorem, we can conclude that:

[0012]

[0013] in, is the average value of blood flow perfused through the VOI. Referring to myocardial perfusion and cardiac output (CO), we obtain the following relationship:

[0014]

[0015] It describes the volume of blood (eg, in liters) pumped by the heart per unit time (eg, per minute) as the amount of contrast agent at the inlet divided by the area under the concentration curve (AUC) recorded at the outlet.

[0016] If we analyze the in vivo recordings of the arterial input function (AIF) and the time attenuation curve (TAC): the change in contrast volume that occurs within a certain time interval is proportional to the difference between the inflow and outflow. This can be expressed as:

[0017]

[0018] It determines the reference volume V of the tissue t The mass flux balance of the contrast agent in t The concentration c(t) in the region monitored by the arterial input AIF curve is denoted by c a (t), the tissue attenuation TAC curve is represented by c t (t), and the concentration of venous outflow is c v (t). Since relative units are used in perfusion analysis instead of absolute values, it is appropriate to express the mean flow Q accordingly by the mean volume of flowing blood in the tissue mass (MBF - myocardial blood flow):

[0019]

[0020] Among them, ρ t is the tissue density. Therefore, the output balance equation takes the form:

[0021]

[0022] When analyzing the shape of the TAC, it can be noted that in the initial phase, the contrast agent gradually accumulates in the control volume, the efflux is negligible, and the concentration curve has a strictly increasing character. For this phase, it is assumed that the analysis is limited to the part with negligible efflux, that is, c v (t)≈0, then the equilibrium equation is in the form of:

[0023]

[0024] And this assumption can be further strengthened by restricting our considerations only to the maximum value:

[0025]

[0026] Therefore, if we want to assess perfusion (MBF - myocardial blood flow) in a practical way, then we can use the maximum slope method. MBF is then expressed as the ratio of the maximum slope of TAC during the accumulation phase to the maximum value recorded on AIF, multiplied by the tissue density ρ t The reciprocal of:

[0027]

[0028] By the same reasoning and by starting from a mass balance instead of its flux in equation (17), we obtain a relation for the unit volume of the vascular bed, which is called myocardial blood volume (MBV).

[0029]

[0030] Calculation of perfusion time parameters (eg, mean MTT and maximum TTP) is a purely arithmetic operation.

[0031]

[0032] TTP is just a coordinate on the timeline, which satisfies the following relationship:

[0033] TAC(TTP)=max(TAC). (25)

[0034] The maximum slope method is an extremely efficient method with moderate hardware requirements. It can be effectively applied even to very large data sets. This is particularly important for analyzing computed tomography (CT) perfusion imaging data.

[0035] The field of non-invasive clinical cardiac imaging is developing rapidly. Developments are focused on quantitative perfusion assessment techniques that can be used to assess coronary blood flow and myocardial ischemia. The main methods for myocardial perfusion assessment are single photon emission computed tomography (SPECT), positron emission tomography (PET), and stress dynamic computed tomography perfusion (CTP) (see AlJaroudi WA, Hage FG 2020 Review of cardiovascular imaging in the Journal of Nuclear Cardiology 2019: Positron emission tomography, computed tomography and magnetic resonance, J Nucl Cardiol, 27(3), pp. 921-930).

[0036] SPECT has been an established method for clinical cardiovascular imaging for decades. Assessment of the extent of myocardial ischemia based on SPECT together with measurements of left ventricular function and volume is considered an accepted method for diagnosing myocardial ischemia. The main goals of SPECT technology development were to maximize sensitivity, reduce acquisition time, and reduce the amount of radiopharmaceuticals injected. SPECT temporal resolution (which is the frame acquisition time) is equal to 10 s, and the voxel size is 10×10×10 mm 3 . In addition, there is also the problem of increased BMI (body mass index) in the patient population. The solution to this problem is cadmium-zinc-telluride solid-state detector technology (CZT), which provides increased sensitivity, increased spatial and energy resolution, reduced acquisition time and smaller radiopharmaceutical doses. CZT SPECT enables dynamic imaging and kinetic assessment of perfusion tracers in the blood and myocardium. Using this information, chamber modeling and depiction of absolute MBF and coronary blood flow reserve are possible (see Dewey M in 2020 Clinical quantitative cardiac imaging for the assessment of myocardial ischaemia, Nature Reviews Cardiology, 17 (7), pp. 427-450).

[0037] Another technique used to confirm or exclude ischemia is PET. PET allows accurate in vivo measurement of radiopharmaceutical concentrations. PET temporal resolution (which is the frame acquisition time) is in the range of 1 s-5 s, and voxel size is 4×4×4 mm 3 . Use different tags in the inspection, but 15 O-water is considered the clinical gold standard.15 Modeling of O-water dynamics is well established, and its tracer concentration is proportional to the measured image signal. Due to the short half-life of the tracer, rest-stress protocols can be completed in 30 minutes, so the total radiation dose is lower than SPECT. Images of perfusion defects are sometimes grainy and difficult to interpret. MBF values ​​can be provided by modeling the tracer dynamics and correcting for limited extraction. Recent advances in PET allow the calculation of parametric images showing voxel-level MBF and facilitate 15 The use of O-water PET in clinical settings (see Dewey M, 2020, Clinical quantitative cardiac imaging for the assessment of myocardial ischaemia, Nature Reviews Cardiology, 17(7), pp. 427-450).

[0038] The last of these three techniques is CTP. There are currently two methods for CTP: static and dynamic. The assessment of myocardial enhancement obtained at a single time point during the first passage of a contrast agent sample is called static CTP imaging. On the other hand, dynamic CTP assesses myocardial perfusion by measuring myocardial enhancement at multiple time intervals during the initial passage of a contrast agent sample, which is inferred to be time-dependent elastic. Dynamic CTP allows a fully quantitative examination of myocardial perfusion. Due to recent advances in computed tomography (CT) technology, dynamic CTP of the entire myocardium can now be performed with significantly less radiation exposure (see Kitagawa K in 2018 Dynamic CT perfusion imaging: state of the art, Cardiovascular Imaging Asia, 2 (2), pages 38-48).

[0039] All the above methods have advantages and disadvantages. PET is superior to SPECT in terms of spatial and temporal resolution, image quality and diagnostic accuracy. PET can also be performed in patients with implanted cardiac defibrillators and pacemakers and in patients with tracers. However, PET involves high costs, high radiation exposure and is not so widely available. SPECT is more suitable for clinical routine because radionuclides are cheaper, easier to prepare and have a longer half-life. In addition, SPECT performs better in detecting ischemia due to its high sensitivity and specificity. SPECT allows the assessment of left ventricular function. The limitations of SPECT are high cost, time consumption, radiation exposure, limited anatomical information due to poor spatial resolution, attenuation artifacts that produce false-positive results (especially in obese patients), and the possibility of underestimating the actual extent of ischemia in patients with multivessel disease. Unlike the other mentioned methods, CTP imaging is the only non-invasive method used to identify the functional significance of coronary artery stenosis and quantify it. A single CTP examination provides a comprehensive anatomical and functional assessment. In addition, it is a very rapid and widely available method with high sensitivity and specificity. CTP allows the detection of smaller perfusion defects due to its high spatial resolution (about 0.3 mm). However, this method suffers from radiation exposure and motion artifacts (respiratory artifacts) from breathing, beam hardening, and high heart rate artifacts (according to Xu C, Yi Y, and Wang Y in Clinical Applications of CT Myocardial Perfusion Imaging, Cardiovascular Imaging Asia, 4 (4) pp. 86-90 in 2020, Seitun S in CT myocardial perfusion imaging: a new frontier in cardiac imaging, BioMed research international, 2018).

[0040] Advances in computed tomography technology in recent years have allowed the clinical application of myocardial perfusion imaging by computed tomography (CTP). In line with PET and magnetic resonance imaging (MRI), CTP can provide qualitative or quantitative perfusion data. Assessment of myocardial perfusion by CT can be achieved through static or dynamic scan acquisition. Although both types of acquisition have been shown to be clinically applicable, important obstacles must be overcome before their routine use. For example, the success of perfusion assessment using static CTP is highly dependent on contrast agent inference timing and motion. In contrast, dynamic CTP for myocardial perfusion assessment allows for a completely quantitative analysis of myocardial perfusion. The myocardium is enhanced at multiple time points during the first pass of the contrast agent, and in this way myocardial perfusion is assessed. Performing stress CTP studies with and without beta-blockers may also affect the assessment of myocardial perfusion (see Seitun S, 2016, Stress computed tomography myocardial perfusion imaging: a new topic in cardiology, Revista de (English edition), 69(2), pp. 188-200).

[0041] Difficulties may arise in establishing clear cutoff MBF values ​​during perfusion studies. Among researchers, cutoff values ​​for CT-derived MBF vary between 75 ml and 164 ml / 100 ml / min. Differences in study design, patient coronary risk profile, sample size, FFR threshold, and CAD prevalence may account for this wide range of cutoff MBF values ​​(see Kitagawa K, 2018, Dynamic CT perfusion imaging: state of the art, Cardiovascular Imaging Asia, 2(2), pp. 38-48). When MBF is uniformly high (i.e., >110 mL / 100 mL / min for dual-source CT and maximum slope techniques), stress myocardial perfusion can be considered normal. Low MBF (75 mL / 100 mL / min) is uncommon but is most likely related to an absent or inadequate hyperemic response to vasodilators, which may make ischemia identification difficult. Once adequate hyperemia is demonstrated by high MBF in the distal myocardium, the relative MBF values ​​in the vascular territory can be used to determine the hemodynamic significance of coronary artery stenosis. Due to normal inter-patient MBF variability, there may not be a generally applicable absolute MBF threshold for myocardial ischemia (see Kitagawa K, 2018, Dynamic CT perfusion imaging: State of the art, Cardiovasc Imaging Asia, 2(2), pp. 38-48).

[0042] With CTP imaging, we can assess what is happening in the body with unprecedented spatial resolution. This helps to better understand myocardial perfusion at the microvascular level. In recent years, developments in this topic have led to a reduction in: scan time, motion artifacts, use of contrast agents, and radiation dose exposure. This has also led to a significant improvement in spatial and temporal resolution. However, all of the above aspects remain problematic and have their limitations.

[0043] For MPI, as mentioned above, we are faced with the problem of artifacts that hinder the accurate assessment of myocardial perfusion, namely artifacts caused by cardiac motion. Artifacts may result in the appearance of hypoperfused areas that may falsely resemble myocardial ischemia. In this case, we are faced with the problem of identifying which defects are related to perfusion and which are false positive findings. The above-mentioned motion artifacts may arise for several reasons, such as due to patient movement during the scanning process or cardiac motion due to high or irregular heart rate. Examining different cardiac phases is a way to reduce the impact of motion artifacts in MPI scans. Specifically, real (or true) perfusion defects can be observed in all phases of the cardiac cycle, while motion artifacts are only visible in a few cardiac phases. In addition, the occurrence of motion artifacts can lead to problems with temporal resolution, especially in patients with high heart rate and high heart rate variability. During dynamic CTP imaging of the heart rate during stress, a sufficiently high temporal resolution as well as a sufficiently high spatial resolution are required to image the entire myocardium "working" at a high rate (according to Seitun S, 2016, Stress computed tomography myocardial perfusion imaging: a new topic incardiology, Revista de (English edition), 69(2), pp. 188-200, Seitun S, 2018, CT myocardial perfusion imaging: a new frontier in cardiac imaging, BioMed research international, vol. 2018, article ID 7295460, pp. 21). Due to their medical conditions, cardiac patients may undergo multiple imaging procedures, which increases their radiation exposure to levels beyond what is expected in a lifetime. For example, dynamic CTP is associated with high radiation doses because this imaging technique requires multiple acquisitions to generate a TAC. During static CTP, the average radiation patient exposure ranged from 1.9-15.7 mSv, with a mean of 5.93 mSv. For dynamic CTP imaging, the effective dose values ​​ranged from 3.8-12.8 mSv, with a mean radiation dose of 9.23 mSv. CT-based dynamic myocardial perfusion imaging requires relatively high radiation exposure. Various techniques are used to reduce radiation dose, for example, setting a low tube voltage during dynamic CTP for patients with a normal body mass index. In this way, the radiation dose can be reduced by 40% while maintaining image quality and good MBF assessment. However, there is still a need to improve the CTP imaging procedure. The improvement is to eliminate unnecessary radiation exposure and the risks associated with the procedure (see Danad I, 2016, Static and dynamic assessment of myocardial perfusion by computed tomography, European Heart Journal-Cardiovascular Imaging, 17(8), pp. 836-844).

[0044] Methods based on indicator dilution flow analysis have been developed for a long time and are not limited to applications in classical medicine. These methods have very wide applications in hydrology and hydrogeology as well as chemical and process engineering. In particular, these methods are used in the case of moving rocks or debris in the riverbed, when the river cross section is uncertain or the water flow is physically inaccessible or unsafe (see Florkowski T, Davis TG in 1969 The measurement of high discharges in turbulent rivers using tritium tracer, Journal of Hydrology, 8 (3), pp. 249-264, Kilpatrick FA, Cobb ED in 1982 Measurement of discharge using tracers, US Department of the Interior, US Geological Survey, Evans GV in 1983 Tracer techniques in hydrology, The International Journal of Applied Radiation and Isotopes, 34 (1), pp. 451-475). People should note the change of terminology between engineering and medicine. In particular, two equivalent engineering terms "dye dilution" and "tracer dilution" were replaced by a medical term "indicator dilution". Thus, cardiac output is measured in physiology and flow is measured in hydraulic engineering. The latter was described by Charles Allen in 1927. Charles Allen named the method the "salt speed" method (see Allen CM in 1927 Hydraulic-turbine tests by the Allen method, Power Plant Engineering, 31 (10), pp. 549-551). The tetrarch of Trachon (37 AD) used chaff to trace underground streams at the possible source of the Jordan River (see Assaad FA, LaMoreaux PE, Hughes T in 2012 Field methods forgeologists and hydrogeologists, Springer-Verlag Berlin and Heidelberg GmbH & Co. KG). From circulatory physiology, the mean transit time (MMT or equivalently tissue transit time - TTT) was introduced.This is because in chemical engineering and reactor design, a very important parameter is residence time or residence time distribution (according to Wingard LB et al. in 1972 Concepts of residence time distribution applied to the indicator-dilution method, J Appl Physiol, 33 (2), pp. 264-275). Nowadays, methods based on indicator dilution flow analysis (also classified as indicator dilution methods) are partially or completely replaced by computational fluid dynamics based on numerical simulation (according to Chen K in 2018 simulation of particle residence time distribution in industrial scale horizontal fluidized bed, Powder Technology, 345 (1), pp. 129-139, Walker P Sheikholeslami R in 2003 Assessment of the effect of velocity and residence time in CaSO4 precipitating flow reaction, Chemical Engineering Science, 58 (16), pp. 3807-3816).

[0045] This review of prior art developments is based on the need for a reliable method that would allow the diagnosis of patients without the above-mentioned risks and problems. The present invention addresses this need.

[0046] Brief description of the invention

[0047] The following description is intended to provide a better understanding of the principles and benefits of the invention as claimed in the appended claims. Therefore, it is not meant to be limiting in any sense. In particular, certain explanations of meaning or function are not to be construed as limiting or exhaustive. They are provided for ease of understanding and cannot replace a complete understanding provided by the general common knowledge of those skilled in the art.

[0048] The present invention redefines the Myocardial Perfusion Imaging (MPI) paradigm. At its core and one of the biggest improvements over the paradigm is the replacement of the in vivo imaging stage (typically achieved with PET, SPECT or MR scanners) with a computer "imaging" stage, while providing medically relevant and consistent results. The central element of the present invention is a model of arterial and myocardial blood flow solved using numerical fluid dynamics techniques. In this way, the shortcomings of all known in vivo methods are eliminated, including all risks to the overall health of the patient undergoing the procedure. In addition to this, a virtually unlimited spatiotemporal resolution is achieved. Furthermore, the results of the present invention are completely free of motion artifacts and the patient is not additionally exposed to radiation. At the same time, the cost of performing such a computer analysis is significantly lower.

[0049] In the scope of the present disclosure, the term “tracker” may be used interchangeably with the term “contrast.” The term “porous” should be understood to include or be equal to the term “saturated porous.”

[0050] In one aspect, the present invention relates to a method for assessing myocardial perfusion in a human patient using non-invasive measurements and permeation mechanics. In particular, the permeation mechanics involves computational fluid dynamics. The method is a computer-implemented method.

[0051] The method includes using patient-specific anatomical data of a human patient to generate an anatomical geometric model of the myocardium connected to a plurality of coronary vessels Patient-specific anatomical data of a human patient is obtained from at least one of the non-invasive measurements. The term "non-invasive measurement" means a measurement that does not involve any type of surgery and / or significant health risk. In some embodiments, the measurement does not include inserting any probe into the patient's body. Generally, the term "non-invasively" means that it does not involve surgery and / or significant health risk. In certain embodiments, "non-invasively" does not include inserting any probe into the patient's body.

[0052] Anatomical geometric model of the myocardium connected to multiple coronary vessels Including fluid domain (Ω f ) and porous domain (Ω p ) and in the fluid domain (Ω f ) and porous domain (Ω p ) interface (Γ). Anatomical geometric model of the myocardium connected to multiple coronary vessels Satisfies the following equation: and This means that the fluid domain (Ω f ) and porous domain (Ω p ) do not overlap. Myocardium connected to multiple coronary vessels is described as and composition. Including fluid domain (Ω f ) and fluid domain (Ω f ) boundary (Γ f ) or by the fluid domain (Ω f ) and fluid domain (Ω f ) boundary (Γ f ), and Contains porous domains (Ω p ) and porous domain (Ω p ) boundary (Γ p )) or by porous domains (Ω p )) and porous domain (Ω p )) of the boundary (Γ p )). In addition, in the fluid domain (Ω f ) and porous domain (Ω p ) satisfies the following equation: This means that the fluid domain (Ω f ) and porous domain (Ω p ) are connected to each other using the interface (Γ). At the same time, the interface (Γ) cannot be "empty". In addition, the interface (Γ) includes porous domains (Ω p ) boundary (Γ p )) and fluid domain (Ω f ) boundary (Γ f ) or by porous domains (Ω p ) boundary (Γ p ) and fluid domain (Ω f ) boundary (Γ f ) composition. Porous domain (Ω p ) can be used to describe blood flow in myocardial tissue. The results of the method allow diagnosis of heart conditions. In particular, the method allows diagnosis of coronary artery disease or heart failure.

[0053] In addition to the above benefits and new methods for performing myocardial perfusion imaging (MPI), the use of fluid domains and porous domains allows for a simplified description of the myocardium connected to multiple coronary vessels. This allows for faster calculations and faster achievement of perfusion metrics. At the same time, it should be emphasized that it may not always be possible to construct an accurate complex model of the myocardium connected to multiple coronary vessels without the use of fluid domains and porous domains.

[0054] In an embodiment, the patient-specific anatomical data is obtained via at least one of computed tomography angiography (CTA), nuclear magnetic resonance (NMR), and / or computed angiography.The anatomical geometric model may be used in computational fluid dynamics to calculate perfusion metrics.

[0055] In an embodiment, the patient-specific anatomical data comprises images stored in DICOM format. This allows reconstruction of the Time Attenuation Curve (TAC) and Arterial Input Function (AIF) curves at multiple time points.

[0056] Generate an anatomical geometry model of the myocardium connected to multiple coronary vessels The step may include segmenting the image included in the patient-specific anatomical structure data to obtain a surface model of the arteries and muscles of the patient's heart. In addition, an anatomical geometric model of the myocardium connected to the plurality of coronary vessels is generated. The step may include the step of discretizing the surface model by an unstructured volume grid.

[0057] The method also includes determining at least one patient-specific boundary condition and selecting at least one physical parameter for describing blood as a flowing medium using patient-specific demographic and general medical data. The boundary conditions and physical parameters will be used in computational fluid dynamics to, for example, calculate perfusion metrics. The term "patient-specific" means data or parameters that are specific to a patient, rather than, for example, for a patient population. The use of patient-specific data allows obtaining results that are medically relevant to a specific patient, rather than average results for a patient population.

[0058] In one embodiment, the patient-specific demographic and general medical data includes at least one of the patient's age, sex, height, weight, systemic blood pressure, blood test results, and / or information about the patient's current medication. The patient's current medication may be selected from, but not limited to, beta-adrenergic blocking agents, angiotensin converting enzyme inhibitors, and / or antiarrhythmic agents.

[0059] In an embodiment, the physical parameter used to describe the blood as a flowing medium is selected from the group consisting of density and dynamic viscosity coefficient.

[0060] The method also includes using computational fluid dynamics to simulate an anatomical geometric model of the myocardium connected to a plurality of coronary vessels The computational fluid dynamics implements the at least one patient-specific boundary condition and the at least one physical parameter for describing the blood as a flowing medium. The implementation allows obtaining patient-specific results useful in diagnosis.

[0061] In a preferred embodiment of the present invention, the fluid domain (Ω f ) can satisfy the following equation:

[0062]

[0063] Where: ρ is fluid density, t is time coordinate, u and p are flow velocity and pressure, and η is the dynamic viscosity coefficient.

[0064] In a preferred embodiment of the present invention, the porous domain (Ω p ) can satisfy the following equation:

[0065]

[0066] Where: ρ-fluid density, u and p-flow velocity and pressure, η-dynamic viscosity coefficient, C F - inertial resistance tensor, K - permeability tensor.

[0067] In a preferred embodiment of the present invention, the fluid domain (Ω f ) and porous domain (Ω p ))'s blood flow can satisfy the following equation:

[0068]

[0069] Among them, for the porous domain (Ω p )) is χ=1, and 0 for other cases.

[0070] For certain embodiments, in the porous domain, the Navier-Stokes equations are supplemented with additional terms associated with the resistance caused by the porous medium. The additional terms are based on the Forchheimer model describing the effects of viscosity and inertial resistance. This combination of the Navier-Stokes equations with the Forchheimer model allows to circumvent the extreme geometric complexity of the porous domain, while it is able to reproduce the nature of the phenomena occurring in said domain. In the method of the present invention, information of key importance for diagnostic purposes, namely information about the flow distribution, is directly extracted. In other words, this unique feature allows to obtain relevant information without the time and effort required to build and solve a mathematical model of the porous domain. This is beneficial in terms of the need for computing power and the speed of calculation.

[0071] The present invention relates to the calculation of a perfusion metric, which in embodiments may include:

[0072] 1) Myocardial blood flow (MBF), which can be defined as the volume of blood flowing through the myocardial tissue at a certain rate,

[0073] 2) myocardial blood volume (MBV), which can be defined as the total volume of blood in a given unit volume of myocardium,

[0074] 3) mean transit time (MTT), which can be defined as the expected time interval that blood spends in a specific part of the myocardium, and

[0075] 4) Time to Peak (TTP), which can be defined as the time interval during which the radiodensity or concentration indicator is at its maximum in a specific part of the myocardium.

[0076] The present invention allows calculation of MBF, MBV, MTT and TTP directly or in the form of indices, ratios, regional values ​​or patient-specific reference values.

[0077] The perfusion metric may be calculated using a regional tracer kinetic recovery method or using sequential Monte Carlo streamline trajectory Lagrangian tracking. In particular, the calculation of the perfusion metric may involve maximum slope (MXS), deconvolution of impulse response function (DRF), multi-pathway, multi-species, and nonlinear blood-tissue exchange (BTX). In other embodiments, the calculation of the perfusion metric may involve:

[0078]

[0079] Among them, s i is the training vector {Δs1,Δs2,…,Δs n}, where Δs i =f(t,r i ).

[0080] Other features and advantages of the present invention will become more apparent from the detailed description and drawings of non-limiting embodiments. It should be understood that, unless otherwise stated, the embodiments presented in the specification and claims can be combined in any order and quantity to produce new embodiments that form part of the present disclosure. This description includes a large number of references to prior art documents, in particular scientific journals. The entire disclosure of the references is part of the present disclosure. BRIEF DESCRIPTION OF THE DRAWINGS

[0082] The present invention will now be described in more detail with reference to the accompanying drawings, in which:

[0083] Figure 1 A process flow diagram incorporating steps of a method according to an embodiment of the present invention is presented,

[0084] Figure 2 Myocardial perfusion calculated using actual tomographic results is presented. Figure 2 Figure (A) shows the tissue attenuation and arterial input curves. Figure 2 (B) shows the long-axis and short-axis views of MBF and TTP.

[0085] Figure 3 is a three-dimensional reconstruction of the cardiac muscle and coronary arteries from computed tomography angiographic projections. Each marker point on the reconstructed three-dimensional grid in the angiographic projections represents an identified stenosis and its counterpart,

[0086] Figure 4 An unstructured tetrahedral 3D mesh of the cardiac muscle and coronary arteries is presented.

[0087] Figure 5 is a visualization of the normalized velocity field (A, B) and pressure field (C, D) calculated using numerical simulation of blood flow, and

[0088] Figure 6 A 17-segment AHA (American Heart Association) model visualization of numerical simulation results of myocardial perfusion for the same human patient (B) compared with reference myocardial blood flow results calculated from CTP (A). The locations of the coronary arteries are marked on both images.

[0089] Detailed Description

[0090] Although this description generally refers to engineering as the source of the medical methods used in the present invention, there are significant differences between engineering methods and medical methods. The differences are visible to those skilled in the art and are not limited to the different nomenclature used in the two fields. For example, certain engineering measurements can be made directly because in a given example, we can measure the flow in a pipe or river, and the shape of the pipe or river does not change significantly during the measurement. On the other hand, the heart is undergoing significant changes in periodicity, but there is no indication at which part of the cycle a given measurement is made. This limited example is only intended to outline the difference between engineering and medicine. This description of the present invention is not meant to explain all aspects of the engineering and medical conventions used herein, as this should be understandable to those skilled in the art.

[0091] The indicator dilution method used in cardiovascular diagnostics, as well as the well-known engineered dye dilution, tracer dilution and even Allen's salt anemometry, all derive from a very basic concept that combines three conventional scales: time, volume and mass. Conventionally labeling them as T, V, M, we have T = V / Q (Q - volume flow rate), M / V = C (C - concentration of the indicator), and given that M varies with time, we have the relationship dM / dt = Q·ΔC (according to Peters AM, 1993, Aunified approach to quantification by kinetic analysis in nuclear medicine, J Nucl Med, 34(4):706-13). The entire analysis with respect to these three conventional scales and the diagnostic parameters (such as MBF, MBV, MTT, TPP) determined by the software attached to the CT scanner is implemented according to equations (1)...(25). The results of in vivo imaging include a series of images, usually saved in DICOM format. The results and their specific format allow the reconstruction of TAC and AIF curves for multiple time points. Each time point is associated with dozens of images that must be saved by the device in order to reconstruct the tissue time decay curve. This may technically require to be done in a split second or even in the time required to record a single cardiac cycle. In fact, the time points on the timeline are quite far away from each other and their mutual time intervals are measured only by tracking the patient's ECG (electrocardiogram). A representative example of this is Figure 2 , Figure 2 The results of medical examinations of a patient are shown: in its upper part (A), the AIF (solid circles) and TAC (open circles) curves and their approximate values ​​for the entire myocardium (here: volume of interest, VOI) can be seen, while in its lower part (B) the distribution of MBF and TTP in the axial plane. Figure 2The maximum slope method implemented as shown in the results (see Equation 22 above) requires knowing max(AIF) and max(dTAC / dt) to determine the MBF. This task is challenging, and due to the large number of perturbations, the first value is determined by approximating the empirical AIF curve (the equation below the horizontal axis). In contrast, finding max(dTAC / dt) is even more challenging. This is usually achieved by determining the angular slope of the tangent line, as shown in the equation in Figure TAC(t) = -102.347 + 21.082t (the tangent line itself is drawn with a thick line). It can be clearly seen that the tangent line is constructed using only three points. It is quite rare to use more than four points for this purpose. As can be seen in the AIF, and even more so in the TAC curve, the in vivo radiodensity curve varies locally from the expected curve. Local variations are not the only variations, because in fact, the entire curve may deviate from the expected curve. It can be seen that the radiodensity values ​​unexpectedly do not decrease, but gradually increase after reaching a local minimum. Although this effect is difficult to observe at a global level, it is a common phenomenon locally for individual voxels and originates from the recycling of the indicator (control) substance.

[0092] The idea of ​​estimating the degree of circulatory insufficiency based on MPI examinations is laudable and will certainly become the new gold standard for diagnostic tests over time. However, its implementation based on developments in engineering (e.g., civil engineering or groundwater hydraulics) does not seem to lead to success. This is explainable. The heart is an organ that is constantly in motion and undergoes periodic contractions. It is impossible to stop the human heart for the needs of diagnostic procedures, even for a few minutes. On the other hand, the technical limitations of diagnostic procedures do not allow a significant increase in temporal resolution. This problem was recognized at a very early stage of the application of perfusion diagnostic methods. For example, a team from the University of Münster concluded: "Sampling intervals greater than 1 second for PCT imaging calculated using software that relies on the maximum slope model significantly changes the absolute CBF and TTP values..." (see Kloska SP et al., 2010, Increasing sampling interval in cerebralperfusion CT: limitation for the maximum slope model, Acad Radiol, 17(1), pp. 61-66). In practice, due to the interrelationship of the cardiac cycle and the time scale of acquisition, and further due to the size of the distribution volume, the derivative of dTAC / dt must be calculated from a segment of the TAC curve, which in turn is constructed using two, three or four time point layers. To make matters worse, the injection site of the contrast agent is usually far from the target organ and the duration of the injection must be sufficiently long relative to the other time scales (as required for the health of the patient).

[0093] In an embodiment of the present invention, it is permitted to use the maximum slope method (MXS), in particular using relations (22)-(25) to recover the regional tracer kinetics (TKR).

[0094] In another preferred embodiment of the present invention, the recovery of regional tracer kinetics (TKR) will be based on the evaluation of how the tracer is retained in the vasculature, which is performed instead of the tracer rate. If Δm(t) of the total amount m of contrast agent introduced at the inlet has flowed into the exit region of interest (ROI) within the time interval Δt, the empirical probability of such an event is expressed by:

[0095]

[0096] In this way, we can construct an empirical probability density function h(t) that relates the probability value P(t∈Δt) to the interval width Δt:

[0097]

[0098] Knowing the distribution of the probability density function h(t) versus the time it takes for a portion of the tracer material to reach the outlet, the empirical cumulative distribution function can be easily determined by using the following equation:

[0099]

[0100] and the probability of the opposite event consisting of an event whose time exceeds an assumed value of t:

[0101] R(t)=1-H(t). (29)

[0102] Mathematically, this can be expressed as follows: if m number of tracers are introduced at the inlet at infinitesimal time intervals, then at any time t, the amount H(t)·m reaches the ROI volume, and there are still R(t)·m tracers left outside the ROI. Under real-world conditions, the tracers are introduced at infinitesimal time intervals Δt. i The fixed part is usually introduced into the distal position and it slowly flows with the blood toward the assumed ROI. Then, by i At the inlet, the earlier injection (which corresponds to the inlet value AIF (t i The concentration change recorded by the TAC(t) curve caused by i )·MBF·ρ t ·Δt i If we consider the gradual change in concentration recorded at the arterial input function as a function of duration Δt1, Δt2, …, Δt nThe time t1, t2, …, t n If a series of injections is given at , the concentration in the ROI volume at time t will be expressed by:

[0103]

[0104] And therefore, with sufficiently large n→∞, the concentration in the target region of interest will become:

[0105]

[0106] By introducing the impulse response function IRF(t) = ρ t MBF R(t), it is easy to see that the resulting product integral defines the convolution operation, which can be expressed as:

[0107]

[0108] The resulting equation is the central element of the secondary and preferred embodiment of the present invention. This equation defines a simple transformation principle of the state observed at the input to the state in any region of interest, and this can be expressed as follows:

[0109]

[0110] Awareness of the impulse response function IRF(t) greatly facilitates the assessment of perfusion parameters (perfusion metrics). Note that since max(R(t)) = 1, then by definition max(IRF(t)) = ρ t MBF, therefore myocardial blood flow MBF = max(IRF(t)) / ρ t By summing equation (31) over a long period of time, we obtain the following equation using Fubini's theorem:

[0111]

[0112] because But MBF = MBV / MTT. The above method is part of DRF.

[0113] In a third preferred embodiment of the present invention, regional tracer kinetics (TKR) are restored by the blood-tissue exchange model (BTX) (according to Wirestam R in 2012, Using contrast agents to obtain maps of regional perfusion and capillary wall permeability, Imaging in Medicine, 4(4), pp. 423-438, Bindschadler M et al. in 2014, Comparison of blood flow models and acquisitions for quantitative myocardial perfusion estimation from dynamic CT, Phys Med Biol, 59(7), pp. 1533-56, Tofts PS et al. in 1999, Estimating kinetic parameters from dynamic contrast-enhanced T(1)-weighted MRI of adiffusable tracer: standardized quantities and symbols, J Magn Reson Imaging, 10(3), pp. 223-32). In the context of pharmacokinetic and pharmacodynamic (PK / PD) modeling, the previously proposed TKR model is a one-compartment model. In fact, it would be easy to construct a corresponding model of the blood-tissue exchange process on top of the n-compartment PK / PD model. For example, a two-compartment exchange model can be defined as follows:

[0114]

[0115] Wherein, c is the tracer concentration (p-plasma, e-extracellular, extravascular space, ap-arterial plasma), v is the fractional plasma volume in a compartment, and PS is the permeability-surface area product. In embodiments where the effect of the tracer in the extravascular space is negligible, c t =v e c e , and dc t / dt=v e dc e / dt.

[0116] The present invention includes computer determination of the perfusion parameters (perfusion metrics) MBF, MBV, MTT and TTP using computational fluid dynamics (CFD) methods. In a specific embodiment of the present invention, the determination is achieved by using one of the regional tracer kinetic recovery (TKR) methods: maximum slope (MXS), deconvolution of impulse response function (DRF), multi-path, multi-species, non-linear blood-tissue exchange (BTX). Assuming the metric TAC(t) as a result of arterial input, the use of any of the listed methods is associated with the familiar assumption of the temporal evolution of the marker concentration in the volume ROI at each point (voxel). In the case of in vivo testing, the temporal evolution of the marker concentration in the volume ROI at each point (or voxel) is directly provided by the imaging results and expressed in Hounsfield units (CFD does not show these same things, CFD provides pressure and velocity distributions). After image co-registration and motion artifact removal, the tissue attenuation curve in the volume of interest (VOI) is reconstructed voxel by voxel. This commonly adopted strategy is used by every manufacturer of perfusion diagnostic equipment because it is the simplest. The simplest strategy cannot cover all relevant situations and conditions. As mentioned previously, factors influencing the measurements are leading, especially in this case where the strategy is usually adopted, to differences between the results obtained with different devices and / or with different protocols and / or with different critical values ​​of MBF, MBV, MTT or TTP. The differences between the results provided by different researchers can be tens of percentage points (sometimes exceeding 100%). In addition, the results obtained using this approach do not have any direct physical interpretation, which can be expressed in terms of, for example, established flow, volume or time scales. In addition, the results are indeed loosely related to the flow, volume or time itself. Therefore, more and more researchers are interested in the topic of in vivo perfusion imaging via advection-diffusion (see Liu P et al., Perfusion imaging: An advectiondiffusion Approach, IEEE Transactions on Medical Imaging, 40(12), pp. 3424-3435, 2021). If any substance is transported at this flow rate, then its local concentration will change, Described as the advection-diffusion transport equation (Bird RB et al., 2006, transport phenomena, John Wiley & Sons), it can be expressed as:

[0117]

[0118] Here, q is the source term, and the flux (more precisely, the flux density is the ratio of quantities per unit area) is given by:

[0119]

[0120] In a preferred embodiment of the present invention, the marker concentration field distribution can be determined based on the above-mentioned advection-diffusion transport equation. More precisely, by numerically solving the equation. In the case of in vivo imaging result analysis, which is contrary to the computer method that is the subject of the present invention, the velocity field does not have to be calculated, but can be relatively easily reconstructed based on the imaging results. It is difficult to estimate how this method will be recognized by perfusion manufacturers as an in vivo imaging device. In any case, this method has a great inspiration for computer methods. In this case, the velocity field cannot be reconstructed, but must be calculated using a separate model. From an engineering point of view, the description of the marker flow in the tissue can be understood as a description of the marker flow in a porous medium or a saturated porous medium. In general, the description of this flow can implement the Henry Darcy equation (see Ingham D, Pop I in 2005, Transport phenomena in porous media, Vafai K in 2015, Handbook of porous media, CRC Press), which has the following form:

[0121]

[0122] Where K is the permeability tensor. Considering that the velocity field satisfies the Euler equilibrium equation (for liquids), the following differential equation can be obtained

[0123]

[0124] This equation gives the possibility to calculate the pressure distribution in the porous domain. Once we reach the pressure (p), we can use Darcy's equation to obtain the velocity u (see Alves JR et al., 2019, Simulation of the perfusion of contrast agent used in cardiac magnetic resonance: a step toward non-invasive cardiac perfusion quantification, Front Physiol, 10, p. 177, Alves JR et al., 2016, Simulation of cardiac perfusion by contrast in the myocardium using a formulation of flow in porous media, J. Comput. Appl. Math, 295 (C), pp. 13-24). Each differential equation case (including equation (36)) requires solving the boundary conditions (here: pressure). In each case, also in a preferred embodiment of the invention, the pressure distribution on the edge of the porous cardiac muscle is affected by the inflow of the coronary arteries. This observation can be extended to perfusion models of different organs. In a simplified approach, the Daniel Bernoulli equation can be used to achieve this goal. (According to Lukes et al., Numerical simulation of liver perfusion: from CT scans to FE model, CoRR, abs / 1412.6412, 2014). A more advanced approach to modeling myocardial inflow was proposed in the work of a French-American team (see Papamanolis L et al., Myocardial perfusion simulation for coronary artery disease: a coupled patient-specific multiscale model, Ann Biomed Eng, 49(5), pp. 1432-1447, Jaquet C, Toward myocardium perfusion from x-ray computed tomography, Ph.D. thesis, University of Paris-Est, 2018), which developed the following equation

[0125]

[0126] The equation describes the coronary blood flow as a stable blood flow along the centerline axis of each blood vessel. In fact, the model that realizes the equation is well known to any engineering research textbook on hydrodynamics (see Wylie EB, Streeter VL in 1967 Hydraulics Transients, McGraw-Hill). It is worth noting that this method of describing coronary blood flow can be easily and effectively expanded by introducing inertial effects. A good example is the work of Siemens Research and Technology Team (see Sharma P et al. A patient-specific reduced-order model for coronary circulation in 2012, 2012 9th IEEE International Symposium on Biomedical Imaging, pp. 832-835, Formaggia L et al. One-dimensional models for blood flow in arteries in 2003, Journal of Engineering Mathematics, 47, pp. 251-276). However, so far, its original form has not been used in the context considered here. Equation 38 is not used in the present invention.

[0127] In the preferred embodiment of the invention, the velocity (u) and pressure (p) fields are described by the Euler mass balance and Navier-Stokes momentum equations. The myocardium (and muscle tissue in general) supported by a network of branching arteries can be considered as two non-overlapping and simply connected domains. in is the fluid domain Ω f The closure of f Fluid domain Ω f ,and Indicates that there is a boundary Γ p The saturated porous medium domain Ω p The interface between the two domains can be represented by Indicates. In the fluid region Ω fIn the paper, blood motion is considered as a transient flow in an incompressible medium described by the Euler mass balance and Navier-Stokes momentum equations of the following form (based on Bird RB et al., Transport phenomena, 2006, John Wiley and Sons, Pozrikidis C, Introduction to theoretical and computational fluid dynamics, Oxford University Press, 2011, Chung TJ, Computational Fluid Dynamics, Cambridge University Press)

[0128]

[0129] As is known, blood forms a polydisperse system consisting of a liquid phase and a suspended dispersed phase. The liquid phase is plasma, and the suspended dispersed phase is produced by the morphogenetic components: red blood cells (which are the most numerous and are responsible for gas exchange), white blood cells (which are less numerous and provide immune functions), and platelets (which actively participate in the coagulation process, ensuring hemostasis). Under average conditions, the density of blood ranges from ρ f =1.05…1.06g / cm 3 Blood density generally decreases with age (reaching a maximum in newborns). Blood density can be determined more precisely as a function of the volume fraction of morphogenetic components present in human blood. For example, Hawksley's formula is often used to estimate blood density: where ρ pl ,ρ RBCIt is the density of plasma and red blood cells, and HCT is the hematocrit (see De Gruttola S et al. in 2005 Computational simulation of anon-newtonian model of the blood separation process, ArtifOrgans, 29 (12), pp. 949-59). Compared with blood density, the estimation of blood viscosity coefficient is a much more complicated task (according to Milnor WR in 1982 Hemodynamics, Williams and Wilkins, Vlachopoulos CH et al. in 2011 McDonald's blood flow in arteries, CRC Press, Kannojiya V et al. in 2021 Simulation of Blood as Fluid: A Review from Rheological Aspects, IEEE Rev Biomed Eng, 14, pp. 327-341). First, the coefficient depends on the content of the solid part, which is considered to be a typical feature of a dispersed system. For blood, the coefficient is expressed as the value of the hematocrit number (HCT). In addition, erythrocytes have a unique ability to generate three-dimensional structures in the form of rouleaux (see Samsel RW, Perelson AS in 1982 Kinetics of rouleau formation, IA mass action approach with geometric features, Biophys J, 37 (2), pp. 493-514, Samsel RW, Perelson AS in 1984 Kinetics of rouleau formation. II. Reversible reactions, Biophys J, 45 (4), pp. 805-2). In a preferred embodiment of the present invention, a Williamson-like model (Williamson fluid model) is used to describe the blood rheological properties to be used in computational fluid dynamics. For this model The structuring function of the F model in the preferred embodiment is selected according to the Carreau, Cross, and Yasuda model (see Carreau PJ, Rheological equations from molecular network theories, Journal of Rheology, 16(1), pp. 99-127, 1972, Cross MM, Rheology of non-Newtonian fluids: A new flowequation for pseudoplastic system, Journal of Colloid Science, 20(5), pp. 417-437, 1965, Yasuda K, Investigation of the analogies between viscometric and linear viscoelastic properties of polystyrene fluids, Ph.D. thesis, Massachusetts Institute of Technology). In another embodiment, the Walburn-Schneck model can also be used, which explicitly takes into account the influence of plasma proteins (see Walburn FJ, Schneck D, A constitutive equation for whole human blood, Biorheology, 1, pp. 201-210, 1976).

[0130] Engineers usually use Darcy's equation (see equation (36)) to express the porous domain Ω p This equation has been explicitly rejected and is not implemented in the present invention. The equation does not include inertial effects in its construction (in principle, it can only be used in the case of Reno p <1). In contrast, in a preferred embodiment of the present invention, a nonlinear equation is used (according to Forchheimer P, Wasserbewegung durch Boden, Z. Ver. Deutsch. Ing, 45 (1901), pp. 1782-1788, Whitaker S, The Forchheimer equation: A theoretical development, Transport in Porous Media, 25, p. 2761, Vafai K, Handbook of porous media, CRC Press, 2015), which has the following form:

[0131] In Ω p Chinese, (40)

[0132] Among them, C F is the inertial resistance tensor. In order to obtain a consistent phenomenological description, in a preferred embodiment of the present invention, a penalty method is introduced to f and porous domain Ω p The flow of the part is modeled (according to BruneauCH, Mortazavi I in 2008, Numerical modelling and passive flow control using porous media Computers & Fluids. 37 (5): pp. 488-498), which has the following equation

[0133]

[0134] When in the porous domain, χ = 1, and otherwise it is 0.

[0135] Since the solution of the above equations provides a complete description about the evolution of the pressure field and about the velocity (in the grid nodes used for calculation) and therefore about the local value of the flow velocity in each grid cell, in another embodiment of the invention the myocardial blood flow (MBF) can be calculated directly.

[0136] In another embodiment of the present invention, the reconstruction of the parameters describing the perfusion can be performed by training the vector {Δs1, Δs2, ..., Δs n This is implemented using a sequential Monte Carlo streamline trajectory Lagrangian tracking of the following form:

[0137]

[0138] Where, Δs i =f(t,r i(See Crane M, Blunt M, Streamline-based simulation of solute transport, Water Resour, Res, 35, pp. 3061-3078, Salamon P et al., A review and numerical assessment of the random walk particle tracking method, J. Contam. Hydrol, 87, pp. 277-305, Friman O et al., Probabilistic 4D blood flow mapping, International Conference on Medical Image Computing and Computer-Assisted Intervention, 13 (Pt 3), pp. 416-23, Bear J, Modeling phenomena of flow and transport in porous media, Springer International, 2018.) This form of post-processing of the calculation results enables the distribution of the probability density function h(t), the distribution function H(t), and the residual function R(t) to be determined according to equations (27)-(29). In this way, any of the above methods can be used to calculate perfusion parameters based on the results of CFD simulation.

[0139] An embodiment of the method of the present invention includes six steps: 1) collection of patient-specific demographic and general medical data as well as patient-specific anatomical data, 2) quantification of boundary conditions, 3) creation of discrete geometric models, 4) numerical simulations in fluid and porous domains, 5) recovery of regional perfusion metrics, and finally 6) visualization of regional perfusion metrics.

[0140] In the first step, patient-specific demographic and general medical data are collected. In particular, patient-specific demographic and general medical data include age, sex, height, weight, systemic blood pressure and / or blood test results (if available, particularly hematocrit number) of a particular patient. The patient's demographic and general medical data may further include information about the patient's medication. The patient's medication may include beta-adrenergic blockers, angiotensin converting enzyme inhibitors and / or antiarrhythmic agents. Other embodiments include other medications. It should be noted that medication also includes different dosing regimens and any further medical effects of medication. It is expected that in certain embodiments, only selected patient-specific demographic and general medical data, such as sex, age or pharmacological therapy, are collected. Patient-specific demographic and general medical data may be collected via conversations with the patient, or may be obtained from a database.

[0141] In addition to the patient's demographic and general medical data for which a physical interpretation is subsequently given, patient-specific anatomical data is also collected. In one embodiment, the patient-specific anatomical data includes the results of a computed tomography angiography (CTA) of the coronary arteries and myocardium. In other embodiments, computed tomography and / or magnetic resonance techniques may be used. The patient-specific anatomical data includes information about the anatomical structure of the heart myocardium and the blood vessels connected to the myocardium. The patient-specific anatomical data may be collected before, after, or in parallel with the collection of the patient-specific demographic and general medical data. The patient-specific anatomical data may be retrieved from a database.

[0142] In a second step, the patient's demographic and general medical data are used to determine one or more patient-specific boundary conditions, as well as to select one or more physical parameters that describe the blood as a flowing medium. This is required for numerical simulations. The selected physical parameters may include density and / or dynamic viscosity coefficient. Figure 1 In the illustration of , this is indicated by 2.A. In a preferred embodiment, at the same time, the imaging results from the CTA (preferably in the form of a series of hundreds of DICOM images) are segmented, which ultimately provides a surface model of the arteries and muscles of the patient's heart. This is preferred because it allows the results of the method to be obtained more quickly. In other embodiments with less computational requirements, the segmentation can be performed before or after determining the one or more patient-specific boundary conditions. Figure 1 In the diagram, this is represented by 2.B.

[0143] As an example, Figure 3The orientation of the main coronary arteries visualized by curved planar reconstruction (CPR) of the arteries along the centerline is shown (CPR view). As is known to the skilled person, curved planar reconstruction is sometimes referred to as curved planar reconstruction. Figure 3 The case shown has two stenoses in the left branch (a: 92% and d: 87%, OM branch) and three stenoses in the right branch (b: 57%, c: 64%, e: 69%). In the upper part of the figure, the CPR view can be seen, and in the lower part of the figure, the result of the segmentation and reconstruction of the surface model is magnified (again, the locations of the stenoses are marked).

[0144] The surface model obtained by Figure 4 The unstructured volume grid shown is discretized. Figure 4In the upper left corner of the figure, we have a complete view of the computational meshes for the fluid and porous domains, and the rest of the figure shows the expanded artery and myocardial meshes. In addition, changes in radiodensity recorded in CTA are used to quantify the perfusion characteristics of the myocardium (according to Vinnakota KC, Bassingthwaighte JB, Myocardial density and composition: a basis for calculating intracellular metabolite concentrations, American Journal of Physiology-Heart and Circulatory Physiology, 286(5), pp. H1742-H1749 in 2004, Gheorghe AG et al., Cardiac left ventricular myocardial tissue density, evaluated by computed tomography and autopsy, BMC Medical Imaging, 19, p. 29 in 2019, Cuong NLQ et al., Porosity estimation from high resolution CT SAN images of rock samples by using Housfield unit, Open Journal of Geology, 8, pp. 1019-1026 in 2018, Sudhyadhom A, On the molecular relationship between Hounsfield Unit (HU), mass density, and electron in 2020 density in computed tomography (CT), PLoS ONE 15(12), page e024486).

[0145] In the third step, the computational grid with conditions set at its edges (boundary conditions) and with fluid specifications and physical properties of the porous domain is transferred to the solver engine. As previously mentioned, two sets of equations are solved, albeit in a penalty method. For the entire computational grid, the fluid motion equations are solved in the form of the Euler and Navier-Stokes equations. Those skilled in the art will appreciate that the entire computational grid excludes the boundary conditions set to achieve a specific solution. In the porous domain, the Navier-Stokes equations are supplemented with additional terms associated with the resistance caused by the porous medium. The additional terms are based on the Forchheimer model describing the viscosity and inertial resistance effects. This combination of the Navier-Stokes equations and the Forchheimer model allows bypassing the extreme geometric complexity of the porous domain, while it is able to reproduce the properties of the phenomena occurring in the domain. In the method of the present invention, information of key importance for diagnostic purposes, i.e., information about flow distribution, is directly extracted. In other words, this unique feature allows relevant information to be obtained without the time and effort required to construct and solve the mathematical model of the porous domain. This is beneficial in terms of computing power requirements and computing speed. At the same time, it should be emphasized that constructing an accurate and complex model of the myocardium connected to multiple coronary vessels may not always be possible without the use of fluid and porous domains.

[0146] In the fourth step, the results obtained using numerical fluid dynamics methods are interpreted diagnostically by algorithms. Note that, regardless of type and complexity, numerical fluid dynamics models only provide information about the distribution of velocity and pressure fields in the flow. These results are not directly useful for medical diagnosis. When assessing perfusion, the general expectation is that the diagnostic method will provide some contractual metric such as MBF, MBV, MTT, TTP or be related to it in some way (such as a relative index value). The fourth step can be performed in three variants, differing in the computational inputs and the results obtained. In the most direct approach ( Figure 1 - Step 4.A), such as Figure 5 As shown, flow and pressure-flow distributions are determined in views and cross-sections of the long axis of the heart. However, this method does not provide information about conventional perfusion metrics such as MBF. Therefore, in a second variant of post-processing of the simulation results, sequential Monte Carlo streamline trajectory Lagrangian tracking ( Figure 1- Step 4.B). While the usual application of this type of post-processing is the visualization of trajectories, here it is used to reconstruct probability distribution functions: density, distribution and residuals. This information can be analyzed using known techniques used to analyze the results of in vivo myocardial perfusion imaging. Although particle tracking is very effective, its accuracy is limited and it does not replace a full advection-diffusion contrast dynamics simulation. Only this method provides complete information about the contrast agent flow in the tissue and is generally equivalent to in vivo methodologies. These results (4.C) will undergo the same post-processing as the particle tracking results (4.B).

[0147] In the fifth step, reconstruction of the tracer kinetics is performed to determine standard perfusion metrics and is correlated with the results obtained in step 4.B or 4.C (results from 4.A are not subjected to this treatment). Any of the tracer kinetics recovery (TKR) strategies discussed in detail above can be used, i.e., maximum slope (MXS), impulse response function deconvolution (DRF), blood-tissue exchange (BTX). The final results are exported for later visualization.

[0148] In the sixth and final step, only the visualization of the previously obtained perfusion metrics is performed. Figure 5 The visualization can be limited to the presentation of normalized values ​​of flow and pressure distribution. If particle tracking or full contrast dynamics simulations are performed, the visualization will be for myocardial blood flow (MBF), myocardial blood volume (MBV), mean transit time (MTT) and time to peak (TTP). Such results are presented in Figure 6As shown in the figure, here additionally they are placed in the system of the American Heart Association 17-segment model (Cerqueira MD, 2002, Standardized myocardial segmentation and nomenclature for tomographic imaging of the heart, A statement for healthcare professionals from the Cardiac Imaging Committee of the Council on Clinical Cardiology of the American Heart Association, Circulation, 105(4), pp. 539-42). Since perfusion indices vary greatly between patients, the index value can be applied together with normalization for inter-individual differences (Kono AK et al., 2014, Relative myocardial blood flow by dynamic computed tomographic perfusion imaging predicts hemodynamic significance of coronary stenosis better than absolute blood flow, Invest Radiol, 49(12), pp. 801-7, Coenen A et al., 2017, Integrating CT myocardial perfusion and CT-FFR in the work-up of coronary artery disease, JACC Cardiovasc Imaging, 10(7), pp. 760-770).

[0149] Any computational step or sub-step of the method according to the present invention can be implemented using a computer or computer program. In some specific embodiments, some or all of the computations are performed using a computer program stored on a computer or on any type of memory device or both. In another embodiment, some or all of the computations for the purpose of the method can be performed remotely, for example using a cloud-based infrastructure, which can include the use of the Internet or a local network.

[0150] Although all the measurement methods and algorithms described herein are intended to define the parameters of the present invention and provide tangible results compared with clinical trials, they are by no means restrictive, but exemplary. Words such as "including" or "comprising" are not restrictive, and for example, when element A includes another element B, element A may include one or more other elements in addition to element B. The use of singular or plural forms does not limit the scope of the present disclosure, and for example, a part of the description indicating that element A includes element B also discloses that multiple element Bs are included in one element A and multiple element A is included in one element B and wherein multiple element A includes multiple element B embodiments. After reviewing the contents of the present disclosure, many other embodiments will be obvious and clear to those skilled in the art. Therefore, the scope of the present invention should be determined with reference to the attached claims together with the entire scope of equivalents to which the claims are entitled.

Claims

1. A method for evaluating myocardial perfusion in a human patient using non-invasive measurements and permeability mechanics, wherein the method comprises using computational fluid dynamics to simulate an anatomical geometric model of the myocardium connected to a plurality of coronary vessels blood flow in a patient, wherein the computational fluid dynamics comprises at least one patient-specific boundary condition and at least one physical parameter for describing blood as a flowing medium, and in, The anatomical geometric model of the myocardium connected to a plurality of coronary vessels Including fluid domain (Ω f ) and porous domain (Ω p ) and in the fluid domain (Ω f ) and the porous domain (Ω p ) between the interface (Γ), wherein the anatomical geometric model of the myocardium Satisfy the following formula: and and in, Including the fluid domain (Ω f ) and the fluid domain (Ω f ) boundary (Γ f ) or by the fluid domain (Ω f ) and the fluid domain (Ω f ) boundary (Γ f ), and The porous domain (Ω p )) and the porous domain (Ω p ) boundary (Γ p ) or by the porous domain (Ω p ) and the porous domain (Ω p ) boundary (Γ p ), and Wherein, the interface (Γ) satisfies the following formula: and wherein the at least one patient-specific boundary condition is determined using patient-specific demographics and general medical data, and the at least one physical parameter for describing blood as a flowing medium is selected using patient-specific demographics and general medical data, and The anatomical geometric model of the myocardium connected to a plurality of coronary vessels Generated using patient-specific anatomical data of a human patient, and wherein the patient-specific anatomical data of the human patient is obtained from at least one of the non-invasive measurements.

2. The method according to claim 1, wherein: The patient-specific anatomical data is obtained via at least one of computed tomography angiography, magnetic resonance imaging, and / or computer angiography.

3. The method according to claim 2, wherein: The patient-specific anatomical data includes images stored in a DICOM format.

4. The method according to any one of claims 1 to 3, wherein: The anatomical geometric model of the myocardium connected to a plurality of coronary vessels The generating includes segmenting an image included in the patient-specific anatomical data to obtain a surface model of arteries and muscles of the human patient's heart.

5. The method according to claim 4, wherein: The anatomical geometric model of the myocardium connected to a plurality of coronary vessels The generation also includes the step of discretizing the surface model by an unstructured volume grid.

6. The method according to any one of claims 1 to 5, wherein: The patient-specific demographic and general medical data include the patient's age, sex, height, weight, systemic blood pressure, blood test results, and / or information regarding the patient's current medication therapy, wherein the patient's current medication therapy is selected from, but not limited to, beta-adrenergic blockers, angiotensin converting enzyme inhibitors, and antiarrhythmic agents.

7. The method according to any one of claims 1 to 6, wherein: The physical parameter is selected from the group consisting of density and dynamic viscosity coefficient.

8. The method according to any one of claims 1 to 7, wherein: The porous domain (Ω p ) is configured to describe blood flow in myocardial tissue.

9. The method according to any one of claims 1 to 8, wherein: The porous domain (Ω p ) is saturated.

10. The method according to any one of claims 1 to 9, wherein: The fluid domain (Ω f ) satisfies the following formula: And / or wherein the porous domain (Ω p ) satisfies the following formula: and / or wherein, through the fluid domain (Ω f ) and the porous domain (Ω p ) satisfies the following formula: In Ω f ∩Ω p middle, Where: ρ-fluid density; t-time; u and p-flow rate and pressure; η-dynamic viscosity coefficient; in Ω p In χ = 1, and in Ω p Outer χ = 0; C F - inertial resistance tensor; and K - permeability tensor.

11. The method according to any one of claims 1 to 10, wherein: The computational fluid dynamics involves using describes the blood rheology, where η is the dynamic viscosity coefficient and HCT is the hematocrit.

12. The method according to claim 11, wherein: Choose from Carreau, Cross, Yasuda or Walburn-Schneck models.

13. The method according to any one of claims 1 to 12, wherein: The method further comprises using the anatomical geometric model of the myocardium connected to a plurality of coronary vessels The results of the blood flow simulation in are used to calculate perfusion metrics.

14. The method according to claim 13, wherein: The calculating of the perfusion metric comprises using at least one of a regional tracer kinetic recovery method, wherein the regional tracer kinetic recovery method comprises maximum slope (MXS), deconvolution of impulse response function (DRF), multi-pathway, multi-species and nonlinear blood-tissue exchange (BTX), and / or wherein calculating the perfusion metric comprises using a sequential Monte Carlo streamline trajectory Lagrangian tracking of the form: Among them, s i is the training vector {Δs1,Δs2,…,Δs n }, where Δs i =f(t,r i ).

15. The method according to claim 14, wherein: The regional tracer kinetic recovery method includes the following formula: and / or and / or where: τ - instantaneous time at the observation interval [0, t]; AIF - arterial input function; IRD - impulse response function; c - tracer concentration (p - plasma, e - extracellular, extravascular space, ap - arterial plasma); v - fractional plasma volume in the compartment; PS - permeability-surface area product; t - time; TAC - time decay curve; AIF - arterial input function; IRF - impulse response function; MBF - myocardial blood flow; and MTT - mean transit time.

16. The method according to any one of claims 13 to 15, wherein: The perfusion metrics include: Myocardial blood flow (MBF), which is defined as the volume of blood transported through the myocardial tissue at a certain rate, Myocardial blood volume (MBV), which is the total blood volume in a given unit volume of myocardium, Mean transit time (MTT), which is the expected value of the time interval that blood spends in a specific part of the myocardium, and / or time to peak (TTP), which represents the time interval during which the radioactivity density or concentration of an indicator in a particular part of the myocardium will be at its maximum, and Therein, MBF, MBV, MTT and / or TTP are calculated as absolute or relative values.

17. A computer-readable [storage] medium comprising instructions which, when executed by a computer, cause the computer to perform the steps of the method defined in any one of claims 1 to 16.

Images