WO2014077233A1 - 血管解析装置、医用画像診断装置、及び血管解析方法 - Google Patents

血管解析装置、医用画像診断装置、及び血管解析方法 Download PDF

Info

Publication number
WO2014077233A1
WO2014077233A1 PCT/JP2013/080502 JP2013080502W WO2014077233A1 WO 2014077233 A1 WO2014077233 A1 WO 2014077233A1 JP 2013080502 W JP2013080502 W JP 2013080502W WO 2014077233 A1 WO2014077233 A1 WO 2014077233A1
Authority
WO
WIPO (PCT)
Prior art keywords
blood vessel
index
time
analysis
blood
Prior art date
Application number
PCT/JP2013/080502
Other languages
English (en)
French (fr)
Inventor
廣畑　賢治
淳一郎 大賀
Original Assignee
株式会社 東芝
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by 株式会社 東芝 filed Critical 株式会社 東芝
Priority to CN201380054742.6A priority Critical patent/CN104736061A/zh
Publication of WO2014077233A1 publication Critical patent/WO2014077233A1/ja
Priority to US14/716,433 priority patent/US11357409B2/en

Links

Images

Classifications

    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/02007—Evaluating blood vessel condition, e.g. elasticity, compliance
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T7/00—Image analysis
    • G06T7/0002—Inspection of images, e.g. flaw detection
    • G06T7/0012—Biomedical image inspection
    • G06T7/0014—Biomedical image inspection using an image reference approach
    • G06T7/0016—Biomedical image inspection using an image reference approach involving temporal comparison
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/02—Detecting, measuring or recording pulse, heart rate, blood pressure or blood flow; Combined pulse/heart-rate/blood pressure determination; Evaluating a cardiovascular condition not otherwise provided for, e.g. using combinations of techniques provided for in this group with electrocardiography or electroauscultation; Heart catheters for measuring blood pressure
    • A61B5/026—Measuring blood flow
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B5/00—Measuring for diagnostic purposes; Identification of persons
    • A61B5/05—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves 
    • A61B5/055—Detecting, measuring or recording for diagnosis by means of electric currents or magnetic fields; Measuring using microwaves or radio waves  involving electronic [EMR] or nuclear [NMR] magnetic resonance, e.g. magnetic resonance imaging
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/02—Arrangements for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03—Computed tomography [CT]
    • A61B6/032—Transmission computed tomography [CT]
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/02—Arrangements for diagnosis sequentially in different planes; Stereoscopic radiation diagnosis
    • A61B6/03—Computed tomography [CT]
    • A61B6/037—Emission tomography
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/50—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment specially adapted for specific body parts; specially adapted for specific clinical applications
    • A61B6/504—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment specially adapted for specific body parts; specially adapted for specific clinical applications for diagnosis of blood vessels, e.g. by angiography
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/52—Devices using data or image processing specially adapted for radiation diagnosis
    • A61B6/5211—Devices using data or image processing specially adapted for radiation diagnosis involving processing of medical diagnostic data
    • A61B6/5217—Devices using data or image processing specially adapted for radiation diagnosis involving processing of medical diagnostic data extracting a diagnostic or physiological parameter from medical diagnostic data
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B8/00—Diagnosis using ultrasonic, sonic or infrasonic waves
    • A61B8/06—Measuring blood flow
    • G—PHYSICS
    • G16—INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR SPECIFIC APPLICATION FIELDS
    • G16H—HEALTHCARE INFORMATICS, i.e. INFORMATION AND COMMUNICATION TECHNOLOGY [ICT] SPECIALLY ADAPTED FOR THE HANDLING OR PROCESSING OF MEDICAL OR HEALTHCARE DATA
    • G16H50/00—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics
    • G16H50/30—ICT specially adapted for medical diagnosis, medical simulation or medical data mining; ICT specially adapted for detecting, monitoring or modelling epidemics or pandemics for calculating health indices; for individual health risk assessment
    • A—HUMAN NECESSITIES
    • A61—MEDICAL OR VETERINARY SCIENCE; HYGIENE
    • A61B—DIAGNOSIS; SURGERY; IDENTIFICATION
    • A61B6/00—Apparatus or devices for radiation diagnosis; Apparatus or devices for radiation diagnosis combined with radiation therapy equipment
    • A61B6/48—Diagnostic techniques
    • A61B6/481—Diagnostic techniques involving the use of contrast agents
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10072—Tomographic images
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/10—Image acquisition modality
    • G06T2207/10132—Ultrasound image
    • G06T2207/10136—3D ultrasound image
    • G—PHYSICS
    • G06—COMPUTING; CALCULATING OR COUNTING
    • G06T—IMAGE DATA PROCESSING OR GENERATION, IN GENERAL
    • G06T2207/00—Indexing scheme for image analysis or image enhancement
    • G06T2207/30—Subject of image; Context of image processing
    • G06T2207/30004—Biomedical image processing
    • G06T2207/30101—Blood vessel; Artery; Vein; Vascular
    • G06T2207/30104—Vascular flow; Blood flow; Perfusion

Definitions

  • the present embodiment relates to a blood vessel analysis device, a medical image diagnostic device, and a blood vessel analysis method.
  • Coronary artery stenosis is a serious lesion leading to ischemic heart disease.
  • coronary angiography CAG
  • CAG coronary angiography
  • FFR myocardial blood flow reserve ratio
  • the FFR substantially matches the ratio of the stenotic distal coronary pressure to the proximal stenotic coronary pressure.
  • a pressure sensor provided at the catheter tip is measured. That is, the measurement of FFR requires catheter surgery.
  • stenosis analysis of coronary arteries is possible with cardiac CT, compared to FFR measurement by catheter surgery, minimal invasiveness, reduced patient burden, and medical costs can be saved.
  • cardiac CT only an index based on the size of the plaque region or stenosis region included in the CT image can be measured in a minimally invasive manner. If the pressure difference before and after stenosis can be measured by structural fluid analysis based on the CT image, the effect of stenosis (or plaque) can be quantified.
  • Clinical evaluation of coronary circulation dynamics includes ultrafast CT, cineangiogram, ultrasound, nuclear medicine imaging including SPECT (single photon emission tomography) and PET (positron emission tomography), MRI (nuclear magnetic resonance) Imaging method) has been developed and introduced, which is useful for diagnosis and evaluation of treatment methods.
  • SPECT single photon emission tomography
  • PET positron emission tomography
  • MRI nuclear magnetic resonance Imaging method
  • the fluid analysis alone cannot take into account the influence of mechanical factors such as the pulsation of the heart, and therefore cannot accurately measure the blood flow distribution and internal pressure distribution.
  • a structure-fluid interaction analysis that considers the influence of mechanical factors has been performed on the heart and vascular system captured by images.
  • structure-fluid coupled analysis it is often difficult to correctly set the boundary conditions of blood vessel inlet and outlet and the material model of blood vessels and plaques in fluid analysis of blood (including contrast media) .
  • the influence of the micro blood vessels on the blood flow may not be considered. For this reason, the analysis result of the structure-fluid coupling analysis may not reproduce the actual blood flow and blood vessel deformation.
  • the blood vessel analysis device sets a target region for analysis in a blood vessel region included in the time-series medical image, and a storage unit that stores time-series medical image data regarding the blood vessels of the subject, and performs the analysis
  • a setting unit that sets a latent variable identification region in a target region
  • a calculation unit that performs image processing on the time-series medical image and calculates a time-series shape index and a time-series shape deformation index of the analysis target region
  • a construction unit that tentatively constructs a mechanical model related to structural fluid analysis of the analysis target region based on the time-series shape index, the time-series shape deformation index, and the time-series medical image
  • At least one of the predicted value of the blood vessel shape index and the predicted value of the blood flow rate index based on the constructed dynamic model matches at least one of the observed value of the blood vessel shape index and the observed value of the blood flow rate index.
  • a identifying unit for identifying the latent variables related to the latent variable identifying region to include
  • An object of the embodiment is to provide a blood vessel analysis device, a medical image diagnostic device, and a blood vessel analysis method that can improve accuracy related to structural fluid analysis of blood vessels (including blood) and image analysis / tracking processing.
  • FIG. 1 is a diagram showing a schematic block configuration of a medical image diagnostic apparatus (X-ray computed tomography apparatus) according to the present embodiment.
  • FIG. 2 is a diagram illustrating an example of a dynamic model related to a target region for structural fluid analysis according to the present embodiment.
  • FIG. 3 is a diagram showing a typical flow of the structural fluid analysis process performed under the control of the system control unit of FIG.
  • FIG. 4 is a diagram showing a block configuration of the image processing apparatus of FIG.
  • FIG. 5 is a diagram schematically showing a cross section orthogonal to the core line of the blood vessel.
  • FIG. 6 is a diagram showing temporal changes in the shape of the blood vessel core line used in the image tracking process by the image analysis / tracking process of FIG.
  • FIG. 1 is a diagram showing a schematic block configuration of a medical image diagnostic apparatus (X-ray computed tomography apparatus) according to the present embodiment.
  • FIG. 2 is a diagram illustrating an example of a dynamic model related to a target
  • FIG. 7 is a diagram illustrating a calculation example of bending deformation and rotational displacement of the blood vessel core line by the image analysis / tracking process of FIG.
  • FIG. 8 is a diagram for explaining the image tracking process by the image analysis / tracking process of FIG. 4, and shows an example of the tracking process between time t and time t + ⁇ t.
  • FIG. 9 is a diagram showing a cross section orthogonal to the core line of the shape model constructed by the mechanical model construction unit of FIG.
  • FIG. 10 is a diagram for explaining assignment of the forced displacement history to the shape model by the dynamic model construction unit of FIG.
  • FIG. 11 is a diagram illustrating another allocation method of the forced displacement history to the shape model by the dynamic model construction unit of FIG.
  • FIG. 12 is a diagram showing another allocation method of the forced displacement history to the shape model by the dynamic model construction unit in FIG. 4, and is a diagram showing an allocation example when the blood vessel shape deformation index is bending.
  • FIG. 13 is a diagram for explaining the calculation of the posterior distribution and the identification of the average internal pressure regarding the load condition (average pressure in the blood vessel) by the hierarchical Bayes model and the Markov chain Monte Carlo method by the statistical identification unit of FIG. FIG.
  • FIG. 14 is a diagram for explaining the calculation of the posterior distribution and the identification of the material model parameter (equivalent elastic modulus of the blood vessel wall) regarding the material model parameter by the hierarchical Bayes model and the Markov chain Monte Carlo method by the statistical identification unit of FIG. It is.
  • FIG. 15 is a diagram showing a display example of a spatial distribution of internal pressure, which is one of the mechanical indices by the display device of FIG.
  • FIG. 16 is a diagram showing a display example of a spatial distribution of flow velocity values, which is one of the blood flow indexes by the display device of FIG.
  • FIG. 17 is a graph relating to the blood pressure at the left coronary artery origin displayed by the display device of FIG. FIG.
  • FIG. 18 is a graph relating to blood pressure in the vicinity of a branch point between LCX and LDA displayed by the display device of FIG.
  • FIG. 19 is a graph relating to a change in blood pressure related to the direction of the core line displayed by the display device of FIG.
  • FIG. 20 is a diagram illustrating another example of assignment of the forced displacement history by the dynamic model construction unit of FIG.
  • FIG. 21 is a diagram showing another example of assignment of the forced displacement history by the dynamic model construction unit of FIG.
  • FIG. 22 is a diagram showing another example of assignment of the forced displacement history by the dynamic model construction unit of FIG.
  • FIG. 23 is a diagram illustrating a modification of the fiber group according to the present embodiment.
  • FIG. 24 is a diagram illustrating an orthogonal cross section of a dynamic model of a thick cylinder according to the present embodiment.
  • FIG. 25 is an enlarged view of the minute sector element of FIG.
  • the blood vessel analysis device is a computer device for analyzing a structural fluid of a blood vessel region included in a medical image generated by a medical image diagnostic device.
  • the blood vessel analyzing apparatus may be incorporated in a medical image diagnostic apparatus, or may be a computer apparatus such as a workstation separate from the medical image diagnostic apparatus.
  • a medical image diagnostic apparatus such as a workstation separate from the medical image diagnostic apparatus.
  • it is assumed that the blood vessel analysis apparatus according to the present embodiment is incorporated in a medical image diagnostic apparatus for specific description.
  • the medical diagnostic imaging apparatus can be applied to any type of diagnostic imaging apparatus equipped with an imaging mechanism for imaging a subject.
  • a medical image diagnostic apparatus according to the present embodiment for example, an X-ray computed tomography apparatus (X-ray CT apparatus), a magnetic resonance diagnostic apparatus, an ultrasonic diagnostic apparatus, a SPECT apparatus, a PET apparatus, a radiotherapy apparatus, etc. are used as appropriate. Is possible.
  • the medical image diagnostic apparatus according to the present embodiment is an X-ray computed tomography apparatus for specific description.
  • FIG. 1 is a schematic block diagram of a medical image diagnostic apparatus (X-ray computed tomography apparatus) according to this embodiment.
  • the X-ray computed tomography apparatus includes a CT mount 10 and a console 20.
  • the CT gantry 10 images the imaging region of the subject with X-rays in accordance with control from the gantry controller 23 of the console 20.
  • the imaging site is, for example, the heart.
  • the CT mount 10 includes an X-ray tube 11, an X-ray detector 13, and a data acquisition device 15.
  • the X-ray tube 11 and the X-ray detector 13 are mounted on the CT mount 10 so as to be rotatable around the rotation axis Z.
  • the X-ray tube 11 irradiates a subject into which a contrast medium has been injected with X-rays.
  • the X-ray detector 13 detects X-rays generated from the X-ray tube 11 and transmitted through the subject, and generates an electric signal corresponding to the detected X-ray intensity.
  • the data collection device 15 reads an electrical signal from the X-ray detector 13 and converts it into digital data.
  • a set of digital data for each view is called a raw data set.
  • a time-series raw data set relating to a plurality of scan times is transmitted to the console 20 by a non-contact data transmission device (not shown).
  • the console 20 includes a gantry control unit 23, a reconstruction device 25, an image processing device 27, an input device 29, a display device 31, and a storage device 33 with the system control unit 21 as a center.
  • the gantry control unit 23 controls each device in the console 20 in accordance with the scanning conditions set by the user via the input device 29.
  • the reconstruction device 25 generates CT image data related to the subject based on the raw data set. Specifically, first, the reconstruction device 25 generates a projection data set by pre-processing the raw data set. Pre-processing includes logarithmic conversion, nonuniformity correction, calibration correction, and the like. Next, the reconstruction device 25 performs image reconstruction processing on the projection data set to generate CT image data.
  • Image reconstruction algorithms include analytical image reconstruction methods such as the FBP (filtered ⁇ back projection) method, successive approximations such as the ML-EM (maximum likelihood expectation maximization) method and the OS-EM (ordered subset ⁇ expectation maximization) method. Existing algorithms such as image reconstruction can be employed.
  • the reconstruction device 25 generates time-series CT image data based on the time-series projection data set.
  • the CT image includes a pixel region related to a blood vessel contrasted with a contrast agent (hereinafter referred to as a blood vessel region).
  • the CT image may be slice data representing a two-dimensional spatial distribution of CT values, or volume data representing a three-dimensional spatial distribution of CT values.
  • the CT image is volume data.
  • Time-series CT image data is stored in the storage device 33.
  • the image processing device 27 constructs a dynamic model based on time-series CT images and executes structural fluid analysis. Details of the processing of the image processing device 27 will be described later.
  • the input device 29 receives various commands and information input from the user.
  • a keyboard, a mouse, a switch, or the like can be used as the input device 29, a keyboard, a mouse, a switch, or the like.
  • the display device 31 displays various information such as CT images and structural fluid analysis results.
  • a CRT display for example, a liquid crystal display, an organic EL display, a plasma display, or the like can be used as appropriate.
  • the storage device 33 is composed of various storage media such as a hard disk device.
  • the storage device 33 stores various data such as time-series projection data and time-series CT image data.
  • the storage device 33 stores time-series CT image data in a medical image file format conforming to the DICOM (digital imaging and communications in medicine) standard.
  • the storage device 33 may store medical data collected by an external device in association with time-series CT image data in a medical image file.
  • the system control unit 21 includes a central processing unit (CPU), a read-only memory (ROM), and a random access memory (RAM).
  • the system control unit 21 functions as the center of the X-ray computed tomography apparatus.
  • the system control unit 21 executes a blood vessel structure analysis process according to the present embodiment by executing a blood vessel analysis program stored in the ROM or RAM.
  • the system control unit 21, the image processing device 27, the input device 29, the display device 31, and the storage device 33 constitute a blood vessel analysis device 50.
  • the blood vessel analyzing apparatus 50 may be incorporated in a medical image diagnostic apparatus (X-ray computed tomography apparatus), or may be a computer apparatus separate from the medical image diagnostic apparatus.
  • the blood vessel analyzing device 50 receives medical data such as time-series CT images from the medical image diagnostic device or PACS (picture archiving and communication systems) via the network. Collect it.
  • the blood vessel analysis device, the medical image diagnostic device, and the blood vessel analysis method according to the present embodiment can analyze blood vessels in every part of the human body such as a cardiovascular vessel, a carotid artery, and a cerebral artery.
  • the analysis target according to the present embodiment is a blood vessel of the heart.
  • Examples of blood vessels of the heart include coronary arteries and aorta.
  • the coronary artery starts from the coronary artery origin of the aorta, travels on the myocardial surface, and enters the intima side from the epicardium side.
  • the coronary arteries branch into myriad capillaries in the intima of the myocardium. After bifurcation, countless capillaries reintegrate to form the great cardiac vein and connect to the coronary sinus.
  • the coronary vasculature has the characteristic that perfusion must be ensured in the mechanical changes of myocardial contraction and relaxation.
  • a characteristic of coronary blood flow is that it flows more when the perfusion pressure decreases in the left ventricular diastole than in the systole where the internal pressure of the coronary artery origin increases due to mechanical blood flow inhibition by myocardial contraction. Therefore, the normal coronary blood flow velocity waveform is bimodal between the systole and the diastole, and the diastole blood flow is dominant. In hypertrophic cardiomyopathy and aortic stenosis, retrograde waves are observed during systole, and aortic reflux is known to exhibit a specific blood flow waveform depending on the disease.
  • the antegrade waveform in the diastole is closely related to the left ventricular dilation function, particularly the left ventricular relaxation.
  • the peak of the diastolic waveform shifts backward, and the deceleration leg tends to be gentle.
  • coronary blood flow during diastole cannot be sufficiently increased during tachycardia, which may promote myocardial ischemia.
  • Coronary blood flow is generated by applying coronary perfusion pressure equal to the aortic pressure to the left and right coronary arteries branching anatomically from the aortic origin, that is, the pressure at the aortic origin where the coronary artery branches.
  • Coronary blood flow resistance as well as driving pressure, which is aortic pressure is important for determining coronary blood flow.
  • a thick coronary vessel of 140 to 180 ⁇ m or more has about 20% of the resistance to coronary vessels, whereas a microvessel of 100 to 150 ⁇ m or less is said to have a lot of remaining resistance components. Therefore, in the absence of so-called coronary stenosis, the resistance value depends on the tonicity of coronary microvessels.
  • Vascular resistance factors include vascular properties, arteriosclerosis, vascular stenosis, blood viscosity, and mechanical factors.
  • Coronary microvessel tonus is defined by vascular properties, myocardial metabolism (myocardial oxygen consumption), neurohumoral factors, mechanical factors, various vasoactive substances as humoral factors, blood viscosity, cardiac hypertrophy, coronary artery Coronary circulatory disturbance is also affected by various lesions including sclerosis.
  • Coronary blood flow pulsation is affected by the pulsation pattern of coronary blood flow, control of intramyocardial blood flow by myocardial contraction, and the response of intramyocardial blood vessels to mechanical stimulation.
  • the mechanism by which myocardial contraction inhibits blood flow includes an increase in intramyocardial pressure, changes in intramyocardial vascular volume, and compression of intramyocardial blood vessels.
  • the blood flow regulating factor in the myocardial diastole includes diastole coronary artery pressure, diastole extravascular force, heart rate, ratio of diastole to cardiac cycle, and myocardial relaxation.
  • the blood vessel analysis device 50 constructs a dynamic model based on time-series CT images, performs structural fluid analysis on the blood vessels of the heart using the dynamic model, Vessel flow index is calculated with high accuracy. In order to calculate a highly accurate mechanical index or vascular flow index, it is necessary to assign a highly accurate latent variable to the dynamic model.
  • the blood vessel analyzing apparatus 50 performs inverse analysis on the initial dynamic model and statistically identifies latent variables. Thereby, the blood vessel analysis device 50 can determine a highly accurate latent variable.
  • the mechanical index means a mechanical index related to the blood vessel wall.
  • Examples of the mechanical index related to the vascular wall include an index related to the displacement of the vascular wall, an index related to stress and strain generated in the vascular wall, an index related to the distribution of internal pressure applied to the vascular lumen, and a material property indicating the hardness of the blood vessel. Classified into indicators. Examples of the index relating to the material property indicating the hardness of the blood vessel include an average slope of a curve indicating the relationship between the stress and strain of the blood vessel tissue.
  • the blood flow index means a hemodynamic index related to blood flowing through a blood vessel. Examples of the blood vessel flow rate index include blood flow rate, blood flow rate, blood viscosity, and the like.
  • the latent variables are, for example, material model parameters (such as Young's modulus and Poisson's ratio) such as blood vessel material constitutive equation or blood material constitutive equation, load condition parameters such as internal pressure distribution applied to the vascular lumen, and structural analysis.
  • material model parameters such as Young's modulus and Poisson's ratio
  • load condition parameters such as internal pressure distribution applied to the vascular lumen
  • structural analysis at least one of a variation distribution parameter related to the uncertainty of a boundary condition parameter of fluid analysis, a time-series shape index, and a shape deformation index.
  • the variation distribution parameter related to the uncertainty of the time-series shape index and the shape deformation index means that the medical image data includes a variation distribution caused by noise of each CT value and an ambiguity of a boundary threshold value of a living tissue.
  • various uncertainties are expressed as probability distributions.
  • uncertainties include, for example, uncertainties in the spatial coordinates of the boundary coordinates and feature points of blood vessels and blood (such as blood vessel bifurcations and contrast medium dispersion), geometric structure parameters (cross-sections perpendicular to the core line) And the uncertainty of the medical image data itself (CT value, boundary threshold value, etc.).
  • the dynamic model is a numerical model for expressing the behavior of blood vessels and blood.
  • the mechanical model has different types depending on the method of structural fluid analysis.
  • the dynamic model is classified into a continuum dynamic model and a simple dynamic model.
  • the continuum mechanical model is used, for example, in a finite element method (FEM) or a boundary element method.
  • FEM finite element method
  • the simple dynamic model is classified into, for example, a material dynamic model based on material dynamics and a hydrodynamic model based on rheology.
  • the initial dynamic model means a dynamic model to which a sampling set (a set of combinations of parameters) related to the parameters of the latent variable obtained from the probability distribution of the latent variable and the variable range is assigned.
  • FIG. 2 is a diagram illustrating an example of a dynamic model M1 related to a target region for structural fluid analysis (hereinafter referred to as an analysis target region).
  • the dynamic model M1 has an aorta region R1, a right coronary artery region R2, and a left coronary artery region R3. Blood flows from the aorta to the right or left coronary artery.
  • the end of the dynamic model M1 on the aortic root side is set as the blood flow entrance, and the right coronary artery region end and the left coronary artery region end are set as the blood flow outlet.
  • a boundary condition is set for each of the inlet and the outlet.
  • the boundary condition regarding the inlet includes, for example, the flow rate of blood flow at the inlet or the pressure due to blood flow, or the rate of change thereof.
  • the boundary condition regarding the outlet includes, for example, the flow rate of blood flow at the outlet or the pressure due to blood flow, or the rate of change thereof.
  • Deformation of the aorta, right coronary artery, and left coronary artery includes mechanical action on the blood vessel wall caused by blood flow, mechanical action on the blood vessel wall due to heartbeat (external force), load condition at the blood vessel cross section boundary, blood vessel wall Depending on various factors such as the material model, the unstressed state of the vessel, and the geometry of the vessel wall.
  • the mechanical action on the blood vessel wall caused by the blood flow includes, for example, an internal pressure caused by the blood flow and a shearing force caused by the blood flow. Due to the internal pressure caused by the blood flow, the blood vessel is deformed in the radial direction of the blood vessel or in the direction perpendicular to the surface of the blood vessel lumen.
  • the mechanical action on the blood vessel wall due to the pulsation of the heart and the shearing force caused by the blood flow are the aorta region R1, the right coronary artery region as load conditions for deformation of the blood vessel such as expansion, contraction, torsion and bending in the blood vessel core direction.
  • the deformation of the blood vessel such as expansion / contraction, torsion, and bending with respect to the blood vessel core direction is expressed by forced displacement (movement vector or rotational displacement) or temporal change of the load vector.
  • the deformation in the radial direction of the blood vessel or the vertical direction of the lumen surface based on the internal pressure caused by the blood flow is assigned as a temporal change in the pressure distribution to the blood vessel lumen.
  • a displacement constraint condition by forced displacement is assigned to the aorta region R1, the right coronary artery region R2, and the left coronary artery region R3.
  • the degree of deformation of the blood vessel shape depends on the material of the blood vessel wall. Therefore, material models are assigned to the aorta region R1, the right coronary artery region R2, and the left coronary artery region R3. The degree of deformation of the blood vessel shape also depends on the no-stress state of the blood vessel. A residual stress distribution of the blood vessel may be assigned as an initial value of the load condition.
  • the element set consisting of a node set and nodes obtained by discretizing the space of the numerical model for blood vessel analysis, there are areas where the analysis conditions such as the material model, boundary conditions, and load conditions are identified, and areas where they are not identified It is divided into.
  • Displacement constraint conditions of forced displacement history are assigned to the nodes in the analysis condition identification area, and displacement constraint conditions of forced displacement history are assigned only to the nodes on the vascular wall surface (outer surface) in the identification area of the material model.
  • a displacement degree of freedom is secured and no displacement constraint is assigned.
  • transformation freedom degree of structural fluid analysis can be suppressed, and it can analyze stably and efficiently.
  • a dummy dummy set for buffering may be provided on the surface of the blood vessel wall, and a forced displacement may be applied to a node on the surface.
  • Difficult 3 A method for identifying boundary conditions of coronary arteries.
  • Difficult 4 Image analysis and structural fluid analysis using blood vessel shapes with variations due to the uncertainty of medical image data. By overcoming these four difficulties, the blood vessel analysis device 50 achieves improved analysis accuracy compared to conventional blood vessel structure fluid analysis in which latent variables are not identified by inverse analysis.
  • FIG. 3 is a diagram showing a typical flow of the structural fluid analysis process performed under the control of the system control unit 21 according to the present embodiment.
  • FIG. 4 is a diagram illustrating a block configuration of the image processing apparatus 27.
  • the system control unit 21 reads a medical image file to be processed from the storage device 33 and supplies it to the image processing device 27.
  • the medical image file includes, in addition to time-series CT image data, pressure value data related to the subject's blood vessel lumen, blood flow index observed value data, and plaque index.
  • an MRI image or an ultrasonic echo image may be used.
  • the time-series CT image data is data representing a three-dimensional spatial distribution of time-series CT values.
  • the time-series CT images include, for example, about 20 CT images for one heartbeat, that is, about 20 heart phases.
  • the system control unit 21 causes the image processing device 27 to perform area setting processing (step S1).
  • step S ⁇ b> 1 the region setting unit 51 of the image processing device 27 sets an analysis target region for structural fluid analysis in a blood vessel region included in a time-series CT image.
  • the analysis target region is set to an arbitrary part of the blood vessel region related to the coronary artery.
  • the region setting unit 51 sets the analysis target region and the identification target region in the blood vessel region by an instruction from the user via the input device 29 or image processing.
  • FIG. 5 is a diagram schematically showing a cross section perpendicular to the blood vessel core line (hereinafter referred to as a blood vessel cross section).
  • the blood vessel has a tubular blood vessel wall.
  • the central axis of the blood vessel wall is called the core wire.
  • the inside of the vessel wall is called the lumen. Blood flows into the lumen.
  • the boundary between the lumen and the blood vessel wall is called the blood vessel inner wall.
  • a tissue around a blood vessel such as a myocardium is distributed outside the blood vessel wall.
  • the boundary between the blood vessel wall and the tissue surrounding the blood vessel is called the blood vessel outer wall. Plaques may develop inside the vessel wall.
  • FIG. 5 is a diagram schematically showing a cross section perpendicular to the blood vessel core line (hereinafter referred to as a blood vessel cross section).
  • the blood vessel has a tubular blood vessel wall.
  • the central axis of the blood vessel wall is called the core wire.
  • the inside of the vessel wall is called
  • plaques are classified into, for example, calcified plaques, cocoon plaques, and the like. Spider plaques are sometimes referred to as unstable plaques. Wrinkled plaque is soft and has a risk of rupturing the inner wall of the blood vessel and seeping out into the blood vessel as a thrombus. Therefore, it is clinically useful to grasp the characteristics of plaque. Plaque properties and existing areas can be identified by a plaque index included in the medical image file. The plaque index can be relatively discriminated by, for example, the magnitude of the CT value normalized with reference to the CT value of bone. However, it is not easy to analyze the deformation characteristics and hardness of the plaque inside the blood vessel.
  • step S2 the system control unit 21 causes the image processing device 27 to perform image analysis / tracking processing (step S2).
  • step S2 the image analysis / tracking processing unit 53 of the image processing device 27 performs image processing on the time-series CT image to calculate a time-series blood vessel shape index and a time-series blood vessel shape deformation index.
  • the image analysis / tracking processing unit 53 calculates a time-series blood vessel shape index by performing an image analysis process on the time-series CT image, and performs a tracking process on the time-series CT image. A series of blood vessel shape deformation indices is calculated.
  • the image analysis / tracking processing unit 53 extracts a blood vessel region from each CT image in the image analysis processing, and relates to a pixel region relating to a blood vessel lumen (hereinafter referred to as a blood vessel lumen region) and a blood vessel wall.
  • a pixel region (hereinafter referred to as a blood vessel wall region) is specified.
  • the image analysis / tracking processing unit 53 uses a plurality of pixels on a region where a cross section perpendicular to the blood vessel core line or a surface perpendicular to the blood vessel lumen surface intersects the blood vessel lumen, the blood vessel wall, and the plaque region as a blood vessel shape index. The three-dimensional coordinates of are specified.
  • the blood vessel morphology index is not limited to three-dimensional coordinates, and may be various geometric indices.
  • the radius and diameter of the blood vessel lumen at a certain angle in the cross section perpendicular to the core line and the direction vector of 0 °, or the average area and average radius for all angles in the cross section, or a plurality of cross sections perpendicular to the core line direction can be used as the geometric index according to the present embodiment.
  • the image analysis / tracking processing unit 53 sets a plurality of tracking points such as a feature point, a feature shape, a representative point, and a pixel by an instruction from the user via the input device 29 or image processing.
  • the tracking point is set in a blood vessel region, a blood region, a contrast agent region, a proton region, or the like.
  • the image analysis / tracking processing unit 53 sets a tracking point set such as a blood vessel bifurcation or a surface feature shape.
  • the displacement of the node on the surface of the blood vessel wall, inside the blood vessel wall, or inside the blood vessel lumen of the dynamic model A temporal change is calculated by interpolation processing or the like, and the calculated temporal change is assigned as a forced displacement.
  • the image analysis / tracking processing unit 53 defines a node on the blood vessel core line in the dynamic model.
  • the image analysis / tracking processing unit 53 extracts deformation related to expansion / contraction, torsion, and bending with respect to the blood vessel core direction from the temporal change in the displacement of the nodes on the blood vessel wall surface, the blood vessel wall, or the blood vessel lumen of the dynamic model. You may express by giving as a forced displacement of the node in the cross section perpendicular
  • the forced displacement data (forced displacement history) at the nodes at each time in the dynamic model is specified as the blood vessel shape deformation index.
  • FIG. 6 is a diagram showing temporal changes in the shape of the blood vessel core line.
  • a time-series medical image includes 20 CT images per heartbeat. That is, it is assumed that CT images are obtained every 5% from 0% to 95% of the cardiac phase.
  • the core line of the blood vessel region is extracted from each CT image by the image analysis / tracking processing unit 53.
  • the form of the core wire changes as the cardiac phase progresses.
  • FIG. 8 is a diagram illustrating an example of the tracking process between time t and time t + ⁇ t.
  • nodes in the dynamic model from P1 to P10 are set on the blood vessel core line and are mechanically connected to the nodes of the mechanical model of the blood vessel on each cross section. However, it is independent of the nodes of the blood dynamic model. It is assumed that the displacement data of the nodes P1 to P20 on the blood vessel core line is calculated based on the displacement data of the blood vessel tracking point by a process such as interpolation, and the forced displacement is set at each node.
  • a local blood vessel region RA defined by the node P13 and the node P14 is considered.
  • the distance between the node P13 and the node P14 in the core line direction is L and the radius of the blood vessel region is r at time t.
  • FIG. 7 is a diagram showing a calculation example of bending deformation and rotational displacement of the blood vessel core wire.
  • the torsion angle may be calculated from a change in the normal direction vector c of the surface composed of the lines a and b.
  • the image analysis / tracking processing unit 53 performs forced displacement of each node on the core line (movement displacement in the three-dimensional space and the core line direction based on the coordinates of the tracking point and the movement vector. ), And a blood vessel shape deformation index is calculated. For example, the image analysis / tracking processing unit 53 calculates the temporal change in the coordinate difference between two adjacent nodes as the expansion / contraction distance ⁇ L in the core line direction. Further, the image analysis / tracking processing unit 53, for each node on each core line, between the node and another node on the blood vessel region cross section including the node (a node in the blood vessel lumen, the blood vessel wall, or the plaque region).
  • the time change of the distance is calculated as the expansion / contraction distance ⁇ r in the radial direction. Further, the image analysis / tracking processing unit 53 calculates, for each tracking point, a twist angle ⁇ in the core line direction of the node on the core line based on the coordinates and movement vectors of a plurality of tracking points in the vicinity of the tracking point. To do.
  • the image analysis / tracking processing unit 53 may calculate the flow velocity or the average flow velocity or the average flow rate in the cross section in the core line direction as the blood flow index by the image tracking of the contrast agent or proton in the blood region.
  • the blood vessel shape deformation index is used as a forced displacement in the dynamic model.
  • the time-series blood vessel shape index is referred to as a shape history
  • the time-series blood vessel shape deformation index is referred to as a forced displacement history.
  • step S3 the system control unit 21 causes the image processing device 27 to perform a construction process (step S3).
  • step S3 the dynamic model construction unit 55 of the image processing apparatus 27 performs shape history (time-series blood vessel shape index), forced displacement history (time-series blood vessel shape deformation index), and time-series medical image (CT image, MRI image).
  • shape history time-series blood vessel shape index
  • forced displacement history time-series blood vessel shape deformation index
  • CT image time-series medical image
  • CT image time-series medical image
  • a dynamic model relating to the analysis target region is provisionally constructed based on DICOM data such as an ultrasonic echo image).
  • the dynamic model is a numerical model related to an analysis target region for performing structural fluid analysis.
  • the dynamic model construction unit 55 first constructs a shape model for solving a dynamic model (mathematical model) based on the medical image and the shape history.
  • the shape model schematically represents the geometric structure of the blood vessel region at each time.
  • the shape model is divided into, for example, a plurality of discretized regions.
  • the vertices of each discretized region are called nodes.
  • the dynamic model construction unit 55 may construct a shape model for each time based on the blood vessel region and the blood vessel shape index included in the medical image for each time, or the blood vessel region included in the medical image for a specific time phase.
  • a shape model for each time may be constructed based on the blood vessel morphology index.
  • the time phase in which the blood vessel corresponding to the analysis target region is most contracted is set as the time phase in the no stress state. Assume no stress.
  • FIG. 9 is a diagram showing a cross section perpendicular to the core of the shape model.
  • the shape model has a blood vessel lumen region and a blood vessel wall region from the core line toward the outside.
  • a plaque region may be provided in the blood vessel wall region.
  • a dummy element of the blood vessel peripheral tissue may be provided outside the blood vessel wall region.
  • the dynamic model construction unit 55 obtains a sampling value related to the latent variable parameter obtained from the probability distribution and the variable range of each latent variable (for example, from a set of combinations of each parameter by the Markov chain Monte Carlo method). Set Sampling to the dynamic model. For example, as shown in FIG. 2, the dynamic model construction unit 55 sets a boundary condition identification target region (hereinafter referred to as a boundary condition identification region) related to the entrance to the end of the aortic region R1 on the aortic origin side. A boundary condition identification region relating to the exit is set at the end of the coronary artery region R2 and the end of the left coronary artery region R3.
  • a boundary condition identification target region hereinafter referred to as a boundary condition identification region
  • a boundary condition identification region relating to the exit is set at the end of the coronary artery region R2 and the end of the left coronary artery region R3.
  • the dynamic model construction unit 55 assigns sampling values related to boundary condition parameters obtained from the boundary condition probability distribution and variable range to each boundary condition identification region.
  • the dynamic model construction unit 55 includes a material model identification target region (hereinafter referred to as a material model identification region) and a load condition identification target region (hereinafter referred to as a material model identification region) in the aorta region R1, the right coronary artery region R2, and the left coronary artery region R3.
  • a load condition identification region assigns sampling values related to the material model parameters obtained from the probability distribution and variable range of the material model to each material model identification region, and obtains from each of the load condition identification regions from the probability distribution and variable range of the load condition.
  • the dynamic model construction unit 55 may assign the residual stress when the flow rate is 0 to the analysis target region as the initial value of the load condition.
  • the dynamic model construction unit 55 may set a region for which the geometric structure is to be identified (hereinafter referred to as a geometric structure identification region) in a region where the geometric structure is uncertain.
  • the geometric structure parameter is a variation distribution parameter related to the uncertainty of the geometric structure or a variation distribution parameter inherent in the medical image data. The variation distribution of each CT value and the boundary threshold value of the living tissue Or the like.
  • the dynamic model construction unit 55 may set a material model in the plaque region. Details of the material model will be described later.
  • the dynamic model construction unit 55 assigns the time-series blood vessel shape deformation index calculated in step S2, that is, the forced displacement history, to the shape model.
  • a shape model to which a latent variable and a forced displacement history are assigned is called a dynamic model.
  • FIG. 10 shows a part of the mechanical model of the blood vessel and blood as the shape model M2, and is a diagram for explaining the assignment of the forced displacement history to the nodes of the mechanical model.
  • FIG. 10 shows a part of the shape model M2.
  • FIG. 10 shows the case where the core wire is located within M2, the core wire may be located outside M2.
  • PN1, PN2 a plurality of nodes PN (PN1, PN2) are set in the shape model M2.
  • a node on the core line is referred to as PN1
  • PN2 a node in a blood vessel or blood dynamic model
  • the shape model M2 is set on the dummy element surface, the blood vessel outer wall, the blood vessel inner wall, the plaque region surface, the plaque region inside, or the blood region.
  • the dynamic model construction unit 55 assigns a forced displacement, that is, a blood vessel shape change index to each node PN1 of the shape model M2 for each time.
  • the dynamic model construction unit 55 connects the adjacent nodes PN1 and PN1 on the core line with beam elements (or rigid elements) EB. Further, the dynamic model construction unit 55 connects the node PN1 and the other node PN2 included in the orthogonal cross section passing through the node PN1 by the beam element EB. The dynamic model construction unit 55 assigns a constraint condition related to the shape displacement direction of each blood vessel shape deformation index to the node PN1 and the beam element EB.
  • the forced displacement includes the expansion and contraction of the surface of the blood vessel wall (or dummy element), the twist of the surface of the blood vessel wall (or dummy element), and the blood vessel wall (or dummy element).
  • Surface bending deformation For example, in a region where the material model and the internal pressure of the blood vessel lumen are not identified, not only forced displacement in the core line direction but also time-series expansion / contraction (displacement) of the blood vessel wall in the radial direction is assigned as the forced displacement history.
  • the area around the area is not assigned a forced displacement history.
  • the forced displacement history is assigned only to the part (for example, the surface node of the dummy element).
  • a time-series blood vessel shape deformation index is assigned to the dynamic model construction unit 55, the node PN1, and the beam element EB as a forced displacement history.
  • FIGS. 11 and 12 are diagrams showing another method for assigning the forced displacement history to the shape model.
  • FIG. 11 shows an allocation example when the blood vessel shape deformation index is twisted
  • FIG. 12 shows an allocation example when the blood vessel shape deformation index is bending.
  • the dynamic model construction unit 55 may directly assign a forced displacement history to the surface of the shape model or an internal node.
  • the image analysis / tracking processing unit 53 calculates a blood vessel shape deformation index such as an expansion / contraction amount, a torsion amount, and a bending amount for the feature point.
  • the dynamic model construction unit 55 directly assigns the calculated blood vessel shape deformation index by interpolating (interpolating or extrapolating) the nodes around the feature points.
  • the image processing apparatus 27 performs a reverse analysis using the dynamic model provisionally constructed in step S3, and statistically identifies latent variables set in the dynamic model.
  • the statistical identification process is performed in step S6 described later.
  • Steps S4 and S5 are provided for calculating the blood vessel shape index and the blood flow index used for the statistical identification process, respectively.
  • the system control unit 21 causes the image processing apparatus to perform a blood vessel stress analysis process on the image processing apparatus 27 (step S4).
  • the vascular stress analysis unit 57 of the image processing device 27 performs vascular stress analysis on the current dynamic model to calculate a predicted value of a time-series vascular shape index.
  • the blood vessel shape index may be any of the above-described blood vessel shape fingers.
  • a cross-sectional shape index of the lumen region or a cross-sectional shape index of the blood vessel wall in the blood vessel core direction may be used.
  • the cross-sectional shape index of the lumen region is at least one of the coordinate value of the target pixel of the lumen region and the geometric structure parameters of the lumen region (such as the radius of the lumen region and the diameter of the lumen region).
  • the cross-sectional shape index of the blood vessel wall region specifically includes at least the coordinate value of the target pixel of the blood vessel wall region, the geometric structure parameter of the blood vessel wall region (the radius of the blood vessel wall region, the diameter of the wall region, etc.)
  • the predicted value means a calculated value of a blood vessel shape index calculated by performing vascular stress analysis on the dynamic model.
  • step S 5 the system control unit 21 causes the image processing device 27 to perform blood fluid analysis processing (step S5).
  • step S ⁇ b> 5 the blood fluid analysis unit 59 of the image processing device 27 performs blood fluid analysis on the tentatively constructed dynamic model to calculate a predicted value of the time-series blood flow index.
  • the blood flow index is a blood flow rate or a flow rate.
  • the blood flow index may be a spatial or temporal average value of blood flow or flow velocity.
  • the predicted value means a calculated value of the blood fluid index calculated by performing blood fluid analysis on the dynamic model.
  • Step S6 the system control unit 21 causes the image processing device 27 to perform identification processing (Step S6).
  • the statistical identification unit 61 of the image processing device 27 acquires the blood vessel morphology in which the predicted value of the blood vessel shape index calculated in step S4 and the predicted value of the blood flow index calculated in step S5 are collected in advance.
  • the parameters of the latent variables of the mechanical model are statistically identified so as to match the observed value of the index and the observed value of the blood flow index.
  • the statistical identification unit 61 includes a first statistical identification unit 61-1 and a second statistical identification unit 61-2.
  • the first statistical identification unit 61-1 statistically identifies the parameter of the latent variable so that the predicted value of the blood vessel shape index matches the observed value of the blood vessel shape index.
  • the second statistical identification unit 61-2 statistically identifies the parameter of the latent variable so that the predicted value of the blood flow index matches the observed value of the blood flow index.
  • the first statistical identification unit 61-1 and the second statistical identification unit 61-2 will be described in order.
  • the first statistical identification unit 61-1 sets a data distribution based on the predicted value and the observed value of the blood vessel morphology index calculated in step S4.
  • the data distribution indicates, for example, a multivariate normal distribution function related to an error between the predicted value and the observed value of the blood vessel shape index.
  • the statistical identification processing unit 61-1 calculates a normal distribution function value related to an error between the predicted value and the observed value for each node or each element in the dynamic model, and calculates the product of each normal distribution function value. Is set as the data distribution.
  • the data distribution may be set individually for each time or may be set collectively for a plurality of times.
  • the first statistical identification unit 61-1 assigns a prior distribution (prior probability distribution) to the latent variables of the dynamic model. Specifically, a prior distribution is assigned to each of the parameters related to the uncertainty of the material model, boundary conditions, load conditions, and time-series form and shape deformation indices. For example, a prior distribution relating to a pressure value related to a blood vessel lumen, which is one of parameters of a load condition, is assigned. The range of values that can be taken by the pressure value (assumed range) can be empirically limited in advance. The first statistical identification unit 61-1 calculates a probability distribution of internal pressure values for each discretized region, that is, a prior distribution, by executing a Monte Carlo simulation related to the internal pressure values within the assumed range.
  • the dynamic model construction unit 55 has observed that the pressure distribution in the core line direction is smooth, that the pressure change with the passage of time is smooth, and that there is no backflow of blood flow.
  • a probability distribution expressed mathematically by a multivariate normal distribution function may be set as a prior distribution that the average pressure change gradient in the core line direction is negative.
  • a Monte Carlo simulation regarding the parameters of the load condition can be executed, and a sampling value of the load condition (latent variable) for setting in the dynamic model can be obtained.
  • the first statistical identification unit 61-1 calculates a posterior distribution (a posteriori probability distribution) by performing statistical identification processing on the prior distribution and the data distribution for each latent variable.
  • the first statistical identification unit 61-1 identifies the parameters of each latent variable from the statistical values such as the mode value and average value of the posterior distribution for each latent variable. For example, in the case of the above-described example, the posterior distribution related to the vascular lumen pressure value is calculated, and the identification value of the vascular lumen pressure value is calculated from the posterior distribution.
  • FIG. 13 is a diagram for explaining the calculation of the posterior distribution and the identification of the average internal pressure regarding the load condition (average pressure in the blood vessel) by the hierarchical Bayesian model and the Markov chain Monte Carlo method.
  • a calcified plaque region and a hook-like plaque region are set in a blood vessel extending from the blood vessel starting portion.
  • the calcified plaque region is set as the material model identification region, and the calcified plaque region is set as the material model identification region.
  • the blood pressure in the blood vessel decreases as the blood vessel progresses along the blood vessel core direction from the blood vessel starting portion.
  • a plurality of nodes are set along the blood vessel core line.
  • the posterior distribution of the intraluminal pressure is calculated in the orthogonal section (node section) including each node, and the mode value of the posterior distribution is identified.
  • the blood vessel shape index calculated in step S2 is used.
  • the processing by the second statistical identification unit 61-2 is the same as the processing by the first statistical identification unit 61-1 except that the index used for calculating the data distribution is different. That is, the second statistical identification unit 61-2 first sets a data distribution based on the predicted value and the observed value of the blood flow index calculated in step S5. Next, the second statistical identification unit 61-2 assigns a prior distribution to the latent variables of the dynamic model. For example, a prior distribution relating to a material model parameter relating to blood vessels, a material model parameter relating to blood, and a material model parameter relating to plaque is assigned. Examples of parameters of these material models include material model parameters such as elastic modulus and parameters related to viscosity in blood constitutive equations.
  • the assumed range and probability distribution of the parameters of the material model can be set empirically in advance.
  • the second statistical identification unit 61-2 sets a probability distribution of the parameter of the material model, that is, a prior distribution for each discretized region, and performs a Monte Carlo simulation on the parameter of the material model according to the probability distribution limited to the assumed range. And a sampling value of the material model parameter (latent variable) for setting to the dynamic model can be obtained.
  • the second statistical identification unit 61-2 calculates a posterior distribution by performing statistical identification processing on the prior distribution and the data distribution for each latent variable, and calculates each posterior distribution from the calculated statistical value of the posterior distribution. Identify the parameters of the latent variable. For example, in the case of the above-described example, the posterior distribution relating to the parameter of the material model is calculated, and the identification value of the parameter of the material model is calculated from this posterior distribution.
  • FIG. 14 is a diagram for explaining a posteriori distribution calculation and identification of a material model parameter (equivalent elastic modulus of a blood vessel wall) regarding a material model parameter by a hierarchical Bayes model and a Markov chain Monte Carlo method.
  • the blood vessel model is assumed to be the same as that in FIG.
  • the posterior distribution of the parameter (for example, equivalent elastic modulus) of the material model of the blood vessel wall is calculated, and the mode value of the posterior distribution is identified.
  • the observed value of the blood flow index is assumed to be, for example, a change in the blood flow sent to the aorta, and the observed value of the blood vessel shape index is changed from the time-series CT image by image processing. It can be used as a value (CFA).
  • the flow rate and flow rate may be calculated by calculating the temporal change in the amount of movement of the feature points by image tracking of the contrast agent after injection of the contrast agent into the coronary artery.
  • the concentration change amount of the specific region in the blood vessel core direction or temporal direction of the contrast agent is acquired, and the flow rate and the value of the concentration change divided by the distance in the core direction direction of each region and the temporal change rate of the concentration change are calculated.
  • the flow rate may be calculated.
  • image tracking of protons is used, and in the case of ultrasonic echoes, the flow rate is calculated by contrast echography or the like.
  • the geometrical A structure may be included.
  • the statistical identification unit 61 pre-distributes a probability distribution such as a normal distribution representing a variation within a predetermined range of the coordinate value of each node in the core line direction or a variation within a predetermined range of the diameter of the analysis target region. It is good to set as.
  • a restriction that the shape of the analysis target area is smooth and that the order of the nodes in the core line is unchanged may be set as a prior distribution. According to the probability distribution limited to these assumptions, the Monte Carlo simulation of the geometric structure parameters can be executed, and the sampling value of the uncertainty parameter (latent variable) of the geometric structure for setting to the dynamic model Can be obtained.
  • both the statistical identification process by the first statistical identification unit 61-1 and the statistical identification process by the second statistical identification unit 61-2 need not be performed. That is, in each step S6, either the statistical identification process by the first statistical identification unit 61-1 or the statistical identification process by the second statistical identification unit 61-2 may be performed.
  • the first statistical identification unit 61-1 statistically identifies the parameters of the latent variable so that the predicted value of the blood vessel shape index matches the observed value of the blood vessel shape index
  • the second The statistical identification unit 61-2 statistically identifies the parameters of the latent variable so that the predicted value of the blood flow index matches the observed value of the blood flow index.
  • the statistical identification unit 61 matches the predicted value of the blood vessel shape index and the predicted value of the blood flow index with the observed value of the blood vessel shape index and the observed value of the blood flow index.
  • the parameters of the latent variable may be identified statistically. Further details of the statistical identification process by the statistical identification unit 61 will be described later.
  • step S6 the dynamic model construction unit 55 of the image processing device 27 sets the parameter of the latent variable calculated in step S6 to the dynamic model.
  • the system control unit 21 determines whether or not the identification end condition is satisfied (step S7). If it is determined in step S7 that the identification end condition is not satisfied (step S7: NO), the system control unit 21 repeats steps S4, S5, S6, and S7.
  • the identification end condition is expressed by whether or not an index for determining the end of identification (hereinafter referred to as an identification end index) reaches a specified value.
  • an identification end index for example, a difference value between a predicted value and an observed value of a blood vessel morphology index can be cited.
  • the system control unit 21 determines that the identification end condition is not satisfied when the difference value is larger than the predetermined value, and determines that the identification end condition is satisfied when the difference value is smaller than the predetermined value.
  • the identification end index may be, for example, the number of sampling points of the Monte Carlo method. In this case, when the number of sampling points is smaller than the predetermined value, the system control unit 21 determines that the identification end condition is not satisfied, and when the number of sampling points is larger than the predetermined value, the identification end condition is Determined to be satisfied.
  • the dynamic model construction unit 55 sets the latest dynamic model at that time as the final dynamic model.
  • the mechanical model construction unit 55 associates the observation value of the blood vessel shape deformation index, the load condition parameter in the final mechanical model, and the material model parameter (hereinafter referred to as a related model). Called a model).
  • the related model is stored in the storage device 33.
  • the related model may be stored in association with patient information, examination information, or the like for ease of search. Note that the observed values of the blood vessel shape index and blood flow index, the load condition parameter in the final dynamic model, and the material model parameter do not necessarily have to be associated with each other in the form of the model. Also good.
  • a latent variable may be tentatively identified using a simple dynamic model, and then a latent variable may be accurately identified using a continuum dynamic model. .
  • the parameters of the latent variable can be converged in a short time.
  • a simple mechanical model an internal pressure and an external pressure are expressed by a material dynamic equation of a thick-walled cylinder.
  • the Hagen-Poiseuille flow and the modified Bernoulli equation may be used.
  • FEM structural fluid analysis can be mentioned. Details of the identification method using a simple mechanical model and the identification method using a continuum mechanical model will be described later.
  • step S7 When it is determined in step S7 that the identification completion receipt is satisfied (step S7: YES), the system control unit 21 may cause the image analysis / tracking processing unit 53 to perform correction processing (step S8).
  • step S8 the image analysis / tracking processing unit 53 performs the structural fluid analysis result (predicted value of the mechanical index and the predicted value of the blood fluid index) performed based on the latent variable obtained by the inverse analysis by the statistical identification method.
  • the shape of the blood vessel region included in the time-series medical image may be corrected so as to match the observed value (the observed value of the mechanical index and the observed value of the blood fluid index).
  • the display device 31 displays a diagnosis result based on the corrected time-series medical image.
  • the blood vessel analysis device 50 can display the diagnosis result in consideration of the final dynamic model.
  • the display device 31 may display on the screen a blood vessel location / region where the identification by the inverse analysis and the observation result do not match by the structural fluid analysis.
  • an image of a cardiac phase in which the behavior of a blood vessel is fast is often blurred, and there are portions or regions with large errors in the blood vessel shape observed by image analysis based on a medical image.
  • the medical image data in the cardiac phase where the blood vessel behavior is relatively stable has little noise.
  • the blood vessel shape in the cardiac phase having a large error can be correctly interpolated by using the dynamic model.
  • a correctly interpolated shape can be displayed together with a variation distribution from the original data. Thereby, the stability of the blood vessel shape display can be ensured and the user can recognize the uncertainty of the shape.
  • step S8 the system control unit 21 causes the image processing device 27 to perform vascular stress analysis processing (step S9).
  • step S ⁇ b> 9 the vascular stress analysis unit 57 of the image processing device 27 performs vascular stress analysis on the final mechanical model to calculate a spatial distribution of predicted values of time-series mechanical indices. Specifically, the predicted value of the mechanical index is calculated for each discretized region.
  • step S10 the system control unit 21 causes the image processing device 27 to perform blood fluid analysis processing (step S10).
  • step S10 the blood fluid analysis unit 59 of the image processing device 27 performs blood fluid analysis on the tentatively constructed dynamic model to calculate a spatial distribution of predicted values of the time series blood flow index. Specifically, the predicted value of the blood flow index is calculated for each discretized region.
  • FFR may be calculated as a mechanical index or a blood flow index.
  • Step S11 the system control unit 21 causes the display device 31 to perform display processing (Step S11).
  • step S11 the display device 31 displays the predicted value of the time-series dynamic index calculated in step S9 and the predicted value of the time-series blood flow index calculated in step S10.
  • the display device 31 displays a time-series dynamic index or a time-series blood vessel flow index as a moving image of a time-series dynamic model in a color corresponding to the predicted value.
  • the display device 31 holds a color table indicating the relationship between various predicted values and color values (for example, RGB).
  • the display device 31 specifies a color value corresponding to the predicted value using a color table, and displays a discretization region corresponding to the predicted value with a color corresponding to the specified color value.
  • FIG. 15 is a diagram showing a display example of the spatial distribution of internal pressure, which is one of the mechanical indices.
  • the display device 31 displays each discretized area constituting the dynamic model in a moving image with a color corresponding to the internal pressure value related to the discretized area.
  • the user can grasp the dynamic index that changes with time and space by color.
  • FIG. 16 is a diagram showing a display example of a spatial distribution of flow velocity values, which is one of blood flow indexes.
  • the display device 31 displays each discretized area constituting the dynamic model in a moving image with a color corresponding to the flow velocity value related to the discretized area.
  • the user can grasp the blood flow index that changes with time and space by color.
  • the internal pressure of the stenotic site is smaller than the internal pressure of the non-stenotic site.
  • the internal pressure is designated as the mechanical index
  • the user can determine the presence or absence of stenosis based on a local color difference on the dynamic model.
  • the flow rate of the stenosis site is smaller than that of the non-stenosis site.
  • the flow rate is designated as the blood flow index
  • the user can determine the presence or absence of stenosis based on a local color difference on the dynamic model.
  • the vascular stress analysis unit 57 may calculate the spatial distribution of the hardness value as a mechanical index based on the identification result of the material model parameter of the plaque region.
  • the display device 31 may display the spatial distribution of the hardness value related to the plaque region on the dynamic model.
  • the display device 31 may display an internal pressure distribution, a stress distribution, and a strain distribution around the plaque region. The user can use these displays to estimate the nature of the plaque and the likelihood of failure.
  • the predicted values of the mechanical index and the blood flow index are not limited to the method of expressing with the color of the discretized region of the dynamic model.
  • FIG. 17 is a graph relating to the blood pressure of the left coronary artery origin.
  • the vertical axis of the graph of FIG. 17 is defined by normalized blood pressure, and the horizontal axis is defined by cardiac phase [%].
  • FIG. 18 is a graph relating to blood pressure near the branch point between LCX and LDA.
  • the vertical axis of the graph of FIG. 18 is defined by normalized blood pressure, and the horizontal axis is defined by cardiac phase [%].
  • FIG. 19 is a graph relating to changes in blood pressure with respect to the core line direction.
  • the vertical axis of the graph of FIG. 19 is defined by the blood pressure ratio, and the horizontal axis is defined by the distance [mm] from the aorta.
  • the display device 31 allows the user to easily grasp these values by displaying the predicted values of the mechanical index and the blood flow index in a graph.
  • step S11 the structural fluid analysis process ends.
  • the forced displacement history is set in the core line portion and the outer wall portion of the shape model, but the setting location of the forced displacement history is not limited to this.
  • the forced displacement history may be set in a blood vessel wall region between the core line portion and the outer wall portion.
  • FIG. 20 is a diagram showing another example of assignment of the forced displacement history, and shows a cross section of the shape model.
  • the forced displacement history is assigned only to the node PN2 on the outer wall portion OW of the shape model, and the node PN3 of the blood vessel wall region RV is assigned. It is better not to assign a forced displacement history.
  • FIG. 20A when the boundary condition and the material model are identified, the forced displacement history is assigned only to the node PN2 on the outer wall portion OW of the shape model, and the node PN3 of the blood vessel wall region RV is assigned. It is better not to assign a forced displacement history.
  • the forced displacement history may be assigned to both the node PN2 of the outer wall portion of the shape model and the node PN3 of the blood vessel wall region RV.
  • a forced displacement history is assigned to the node PN1 on the core line.
  • the node PN2 and the node PN1 of the outer wall portion OW may be connected by the beam element EB, and the forced displacement history may be assigned to the nodes PN2 and PN3 on the beam element EB.
  • contraction and expansion in the radial direction are expressed by expansion and contraction displacement of the beam element EB.
  • the forced displacement history need not be assigned to the blood vessel lumen region RI.
  • FIG. 21 is a diagram showing another example of assignment of forced displacement history, and shows a cross section of a shape model including a dummy element RD of a tissue around a blood vessel.
  • the dummy element RD is set outside the blood vessel wall region RN.
  • the node PN4 is set in the dummy element RD in addition to the blood vessel wall region RN.
  • a forced displacement history is also assigned to the node PN4.
  • the mechanical model construction unit 55 assigns the forced displacement history to the node PN3 included in the blood vessel wall region RV, and does not assign the forced displacement history when the boundary condition and the material model are not identified. Also good.
  • the forced displacement history is assigned to the node PN3, the material model is identified with reference to the shape index related to the blood vessel wall region RV in addition to the shape index related to the lumen region RI.
  • FIG. 22 is a diagram showing another example of assignment of the forced displacement history, and shows a cross section of the shape model including the plaque region RP.
  • the plaque region RP is included in the blood vessel wall region RV.
  • the plaque region RP is set in the material model identification region.
  • a material model is identified in consideration of the lumen shape index, the blood vessel wall shape index, and the plaque index.
  • the plaque index is, for example, data relating to the plaque property obtained by the tissue property diagnosis by the ultrasonic diagnostic apparatus.
  • the mechanical model construction unit 55 divides the plaque region into a plurality of partial regions according to the properties, and individually sets the material model identification regions in the plurality of partial regions.
  • a parameter range corresponding to the property of the partial area may be set in advance for each partial area.
  • the material model parameters for each partial region are identified by the statistical identification unit 61.
  • the display device 31 displays an index related to the material characteristics of the blood vessel as a mechanical index, so that the user can accurately and easily grasp the property of the plaque.
  • an elastic model As a blood vessel material model, an elastic model, a superelastic model, an anisotropic hyperelastic model, a superelastic model considering viscosity characteristics, and the like are applicable.
  • an anisotropic hyperelastic model for example, a mathematical model proposed by YCFun or a mathematical model called a Holzapfel-Gasser constitutive equation can be applied.
  • the strain energy per unit reference volume is expressed by the following equation (1).
  • the first term of equation (1) represents the energy related to the shear deformation of the isotropic basic material not containing collagen
  • the second term represents the energy related to the volume deformation of the isotropic basic material not containing collagen
  • Item 3 represents the contribution of each group of collagen fibers (considering dispersion in the fiber direction).
  • FIG. 23 is a diagram illustrating deformation of a fiber group. As shown in FIG. 23, a cylindrical outer membrane is assumed. The deformation of the fiber group in the average direction A on the surface defined by the core line direction z and the circumferential direction ⁇ is expressed by the following equation (2).
  • a Casson constitutive expression such as the following expression (3) and an HB constitutive expression such as the expression (4) are suitable.
  • the statistical identification unit 61 uses a statistical method based on a hierarchical Bayesian model and a Markov chain Monte Carlo method as a statistical identification method of a latent variable in a situation where such uncertainty exists.
  • the statistical identification unit 61 sets a data distribution based on the predicted value and the observed value of the blood vessel shape index or blood flow index calculated in step S4.
  • the data distribution indicates, for example, a multivariate normal distribution function related to an error between a predicted value and an observed value of a blood vessel shape index or a blood flow index.
  • the data distribution may be set individually for each time or may be set collectively for a plurality of times.
  • the statistical identification unit 61 assigns a prior distribution to the forced displacement and the latent variable of the shape model.
  • the prior distribution indicates a probability distribution of possible values.
  • the statistical identification unit 61 executes a parameter survey of a numerical simulation regarding the latent variable, and constructs a model that expresses the relationship between the latent variable and the blood vessel shape index or the blood flow index.
  • a model that expresses the relationship between the latent variable and the blood vessel shape index or the blood flow index.
  • the relationship between the material model parameter, the internal pressure distribution parameter, and the blood vessel shape index or blood flow index is defined in the model.
  • the relationship between the blood vessel shape index or blood flow index and the latent variable may be defined by a database or a table instead of a model. These models, databases, or tables are stored in the storage device 33.
  • the statistical identification unit 61 calculates the probability distribution of the blood vessel shape index or the blood flow index from the prior distribution.
  • the statistical identification unit 61 statistically identifies latent variables from the posterior distribution obtained by the hierarchical Bayesian model and the Markov chain Monte Carlo method.
  • This identification problem specifically becomes a defect setting problem that does not satisfy the following three conditions.
  • the three conditions are (1) the existence of a solution is guaranteed, (2) the solution is uniquely determined, (3) the solution continuously changes with respect to the data, and the solution is stable with respect to the measurement error, It is.
  • the defect setting problem is easy to handle when viewed within the framework of normalization theory and its extensions.
  • the standard normalization theory is insufficient to solve the defect setting problem.
  • a methodology for detecting discontinuities in the internal state and using the detected discontinuous points for estimation of the internal state is required.
  • Markov random field theory is effective in the statistical identification processing according to the present embodiment.
  • the statistical identification unit 61 identifies the probability distribution parameter of the latent variable under an appropriate constraint condition. In order to determine appropriate constraints, it is necessary to know the nature of the solution in advance.
  • the statistical identification unit 61 generates a database related to the constraint condition of the solution space based on the simulation and the observed value.
  • the statistical identification unit 61 uses the generated database to execute statistical identification processing based on Markov random field theory and a hierarchical Bayesian model for a super-multi-degree-of-freedom large-scale problem. In the setting of the prior distribution as a constraint condition, the probability distribution of parameters relating to these factors is individually configured in parallel based on many numerical experimental results.
  • the statistical identification unit 61 identifies parameters of latent variables by integrating a plurality of probability distributions and interpolating data loss. For this processing, the statistical identification unit 61 performs estimation by a hierarchical Bayesian method based on a model using Markov random field theory. This is a mechanism that can estimate pressure and flow rate distribution under certain load conditions and boundary conditions based on the identified intermediate variables from the measurement results of the deformation state of the structure to be analyzed.
  • the identification problem of material models, boundary conditions, and load conditions in structural fluid analysis of coronary arteries is positioned as nonlinear inverse analysis, and the uniqueness and stability of the solution are often not guaranteed. Since realistically possible ranges of material properties and blood pressure of living tissue can be assumed as a priori information, these can be set as a probability distribution of prior distribution. In addition, since pressure and displacement can be assumed to be spatially and temporally smooth, this information can also be set as a priori probability distribution as a priori information. Alternatively, when the fact that no back flow occurs in the blood flow can be taken into account, the fact that the gradient of the overall pressure distribution in the blood vessel core direction is negative (there is a pressure drop) may be used as a constraint condition.
  • Square error distribution between observed value of blood vessel shape deformation index based on time-series CT image and predicted value of blood vessel shape deformation index based on dynamic model for load condition (internal pressure distribution, etc.), boundary condition, and material model Can be set as the data distribution.
  • a square error distribution related to the observable average flow rate may be added as a data distribution.
  • a posteriori distribution can be calculated using the hierarchical Bayes model and the Monte Carlo method based on these prior distributions and data distributions.
  • the identification value of the parameter of the latent variable can be obtained by the occurrence probability and variance of the posterior distribution. It can be said that the higher the probability of occurrence and the smaller the variance, the higher the confidence value.
  • an identification value having a small variance among a plurality of identification values may be selected.
  • structural fluid analysis can be performed under each identification condition, each possibility can be recognized, and the identification values and analysis results can be used as guide information for diagnosis and prevention .
  • an error is included in the time-series CT image, an error is also included in the blood vessel shape index related to each node of the dynamic model. For this reason, for example, each blood vessel shape index is treated as a normal distribution random variable having an average value of the predicted value of the blood vessel shape index measured from a time-series CT image, and includes a restriction that the spatial order of positions is maintained.
  • a prior distribution may be set.
  • the dynamic model construction unit 55 can construct different types of dynamic models depending on the type of the dynamic model.
  • the mechanical model construction unit 55 includes both a shape model for stress analysis of the blood vessel wall (FEM model) and a shape model for fluid analysis of blood (FEM model). Build up.
  • the dynamic model construction unit 55 identifies the blood vessel lumen shape, the blood vessel wall surface shape, and the cross-sectional center on the cross section passing through the nodes discretely arranged on the core line.
  • the dynamic model construction unit 55 calculates an average area, an average radius of the lumen, and an average wall thickness based on the blood vessel lumen shape and the blood vessel wall surface shape. Based on the average area, the average radius of the lumen, and the average wall thickness, the dynamic model construction unit 55 constructs a shape model by performing thick-walled cylinder approximation on the blood vessel region of each discretization region.
  • the modified Bernoulli equation in fluid mechanics or Hagen-Poiseuille Flow (Hagen-Poiseuille flow) is used to approximate the average pressure and flow rate of blood flow.
  • a shape model for approximating the relationship between the blood pressure difference and the flow rate of each of the plurality of discretized regions is constructed.
  • the blood pressure difference is a pressure difference between the inlet pressure and the outlet pressure
  • the flow rate means the inlet flow rate (or flow velocity) and the outlet flow rate (or flow velocity) per unit time.
  • the dynamic model construction unit 55 associates the movement vector for each node with the rotational displacement for the node adjacent to the node at each time based on the expansion / contraction distance and the twist angle in the core line direction calculated in step S2. And good.
  • the kinematics of blood vessels and blood are independent of the forces that cause that motion.
  • the basic concept of kinematics for blood vessels and blood is an abstraction of mathematical concepts into intuitive concepts about position, time, objects, motion, and a collection of deformable substances.
  • the basic kinematic tensors that govern the local analysis of deformation and motion for blood vessels and blood are the deformation gradient tensor F and the velocity gradient tensor L.
  • the deformation gradient tensor F defines the size and shape changes that occur in moving blood vessels and blood material elements.
  • the deformation gradient tensor F is represented by the product of a rotation tensor (strict orthogonal tensor) R and a stretching tensor (positive symmetric tensor) U and V.
  • the stretching tensors U, V are brought about by first imposing a stretch on the direction defined by the orthonormal vector R in the basic form and then imposing a rigid body rotation given by the orthonormal vector R. Note that the order of imposing stretch and the order of imposing rigid body rotation may be reversed.
  • the velocity gradient tensor L does not depend on the reference form, but only on the current form.
  • the velocity gradient tensor L defines the rate at which the size and shape changes that occur in the moving blood vessels and blood material elements occur.
  • the velocity gradient tensor L can be separated into a strain rate tensor D (symmetric tensor) and a spin tensor (antisymmetric tensor).
  • the strain rate tensor D represents the rate of change of the stretch as the object just passes through its current form.
  • the spin tensor represents the rate of change of rotation when an object just passes through its current form.
  • Displacement constraint conditions are assigned to deformation gradient tensors and velocity gradient tensors for a part of the outer surface nodes (or integration points) of the mechanical model of the blood vessel.
  • Predicted values deformation gradient tensor, velocity gradient tensor, or their function values (for example, displacement or area)) and observation values (deformation obtained from observation data) regarding the lumen node (or integration point) of the dynamic model
  • the material model parameters, loading conditions (surface force vector in the lumen of the dynamic model), and boundary conditions (force vector at the blood vessel boundary) are matched so that the gradient tensor, velocity gradient tensor, or their function values match. Identify.
  • the initial state of the stress inside the blood vessel may be assumed in advance or may be identified.
  • a dynamic model based on continuum mechanics is based on equations representing the balance of mass, kinematics, momentum, angular momentum, and energy in moving blood vessels and blood.
  • the concepts of mass, force, heat, and internal energy are fundamental.
  • the balance law is that the momentary total derivative of momentum is equal to the sum of object force and contact force, the momentary total derivative of angular momentum is equal to the sum of object torque and contact torque, and kinetic energy and internal It means that the change with time in energy is equal to the sum of power (mechanical energy), heat supply per unit time and heat flux.
  • field equations such as deformation gradient tensors, velocity gradient tensors, and stress tensors can be derived to describe mechanical models of blood vessels and blood.
  • the strain fields in blood vessels and blood satisfy the conformance conditions.
  • the constitutive equation gives a relationship of a set of 10 scalar equations including density, internal energy, velocity vector, stress tensor, heat flux vector, and temperature.
  • the field balance equation gives 8 scalar relationships and added the temperature Give the remaining unknown relationships.
  • the object force b and the heat source r are known.
  • the parameters of the equations giving these scalar fields are material model parameters.
  • the mechanics model based on continuum mechanics is based on a numerical analysis method based on the finite element method or boundary element method, and approximates the field such as displacement vector, stress tensor, strain tensor, and velocity vector, given boundary conditions, load conditions, And can be calculated based on the material model.
  • the method of solving the structure and fluid equations may be either a monolithic method or a partitioned method. Further, the coupling at the interface between the structure and the fluid may be weakly coupled or strongly coupled. In fluid analysis, moving boundaries may be handled by boundary surface tracking methods such as the ALE method, and boundary surface supplemental types such as Immersed Boundary method, Immersed Finite Element Method, Fictitious Domain Method, etc. A technique may be used.
  • FIG. 24 is a diagram showing an orthogonal cross section of a dynamic model of a thick cylinder
  • FIG. 25 is an enlarged view of the minute sector element of FIG. Inner radius r a, the stress and strain when the outer radius r b of the thick cylinder internal pressure p a and the external pressure p b is applied, illustrating the expression of such displacement.
  • E and ⁇ are material model parameters.
  • E represents the elastic modulus and ⁇ represents the Poisson's ratio.
  • the shear stress ⁇ r in the r ⁇ plane can be expressed as the following equation (18).
  • the velocity distribution u is represented by an axisymmetric rotating paraboloid as shown in the following equation (21).
  • the flow rate Q can be expressed as the following equation (22) by integrating the velocity distribution u over the entire pipe cross section.
  • the cross-sectional average flow velocity v can be expressed as the following equation (23).
  • Equation (24) means that the flow rate Q is proportional to the pressure loss ⁇ p.
  • Equation (25) means that the flow rate Q is proportional to the pressure loss ⁇ p. The flow that satisfies this relationship is called the Hagen-Poiseuille flow.
  • Loss heads include friction loss heads on the inner wall of the blood vessel, loss heads such as constrictions (loss heads for diverging and narrowing pipes), branching and bending pipes, and inlet / outlet losses at the inlet / outlet of the blood vessel analysis target area.
  • loss heads such as constrictions (loss heads for diverging and narrowing pipes), branching and bending pipes, and inlet / outlet losses at the inlet / outlet of the blood vessel analysis target area.
  • the pipe friction loss head is expressed by the following equation.
  • ⁇ is the tube friction coefficient.
  • the pipe friction coefficient ⁇ is determined by the Reynolds number Re when the flow is laminar and by the Reynolds number Re and the surface roughness when the flow is turbulent. Frictional resistance is always affected by the viscosity of the fluid. This frictional resistance is energy loss because it consumes power or energy to drive the flow.
  • the image processing device 27 can identify a latent variable using the material dynamics equation, the Hagen-Poiseuille flow, and the modified Bernoulli equation.
  • a dynamic model consider a case in which the deformation of a blood vessel is expressed by an equation of material dynamics of a thick-walled tube, and the change in the tube diameter is expressed by an internal pressure change and an elastic modulus.
  • the stress-free state is the initial shape (for example, the state in which the blood vessel contracts the most)
  • the elastic modulus of the blood vessel wall and plaque is set to a certain value
  • the observed value of the blood vessel shape deformation index such as the average radius of the blood vessel lumen
  • the observed value of the blood vessel shape deformation index is measured from a time-series CT image.
  • the temporal change of the internal pressure distribution of the blood vessel is determined so as to match the temporal change amount of the observed value of the blood vessel shape deformation index.
  • the predicted value of the blood vessel flow rate index is measured by performing blood fluid analysis under this internal pressure distribution.
  • the image processing device 27 changes the initially determined elastic modulus of the blood vessel wall or plaque, and performs the same analysis. By repeating this, the image processing device 27 can detect potentials such as the elastic modulus of the blood vessel wall and plaque that match the observed value of the blood vessel shape deformation index and the observed value of the blood flow index, the internal pressure distribution, and the pressure boundary condition of the fluid analysis. Variables can be determined. In order to perform this determination method more efficiently and stably, a statistical identification method using a hierarchical Bayesian model and a Markov chain Monte Carlo method may be used.
  • the blood vessel analysis device 50 includes the storage device 33, the region setting unit 51, the image analysis / tracking processing unit 53, the dynamic model construction unit 55, and the statistical identification unit 61. .
  • the storage device 33 stores time-series medical image data related to the blood vessels of the subject.
  • the region setting unit 51 sets an analysis target region in a blood vessel region included in a time-series medical image.
  • the image analysis / tracking processing unit 53 performs image processing on the time-series medical image and calculates a time-series shape index and a time-series shape deformation index of the analysis target region.
  • the mechanical model construction unit 55 provisionally constructs a mechanical model related to the structural fluid analysis of the analysis target region based on the time-series form index, the time-series shape deformation index, and the time-series medical image.
  • the statistical identification unit 61 analyzes the region to be analyzed so that the predicted values of the blood vessel shape index and the blood flow rate index based on the tentatively constructed dynamic model match the pre-measured values of the blood vessel shape index and the blood flow rate index. To identify the latent variables of the dynamic model.
  • the blood vessel analysis device 50 identifies latent variables such as a material model, boundary conditions, load conditions, and geometric structures by inverse analysis using a blood vessel shape deformation index and a blood flow index. can do.
  • the blood vessel analysis device 50 repeatedly performs the inverse analysis while changing the latent variable, thereby achieving the above four points of difficulty, namely, 1. 1. Identification method of coronary artery material model 2. Incorporating the effects of heart shape deformation on coronary arteries. 3. Identification method of coronary artery boundary condition; It is possible to identify a latent variable in consideration of all of a material model using a blood vessel shape having uncertainty, a load condition, and a boundary condition identification method. Therefore, the blood vessel analysis device 50 can perform structural fluid analysis in consideration of the influence of external factors such as blood vessels and the heart that are not depicted in the CT image.
  • the accuracy of the structural fluid analysis of blood vessels can be improved.

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Medical Informatics (AREA)
  • Physics & Mathematics (AREA)
  • General Health & Medical Sciences (AREA)
  • Public Health (AREA)
  • Biomedical Technology (AREA)
  • Pathology (AREA)
  • Nuclear Medicine, Radiotherapy & Molecular Imaging (AREA)
  • Heart & Thoracic Surgery (AREA)
  • Veterinary Medicine (AREA)
  • Biophysics (AREA)
  • Molecular Biology (AREA)
  • Surgery (AREA)
  • Animal Behavior & Ethology (AREA)
  • Radiology & Medical Imaging (AREA)
  • High Energy & Nuclear Physics (AREA)
  • Optics & Photonics (AREA)
  • Physiology (AREA)
  • Theoretical Computer Science (AREA)
  • Computer Vision & Pattern Recognition (AREA)
  • Vascular Medicine (AREA)
  • Hematology (AREA)
  • Cardiology (AREA)
  • Pulmonology (AREA)
  • Dentistry (AREA)
  • Oral & Maxillofacial Surgery (AREA)
  • Quality & Reliability (AREA)
  • General Physics & Mathematics (AREA)
  • Data Mining & Analysis (AREA)
  • Primary Health Care (AREA)
  • Databases & Information Systems (AREA)
  • Epidemiology (AREA)
  • Apparatus For Radiation Diagnosis (AREA)
  • Magnetic Resonance Imaging Apparatus (AREA)
  • Image Processing (AREA)
  • Image Analysis (AREA)

Abstract

　血管の構造流体解析の精度の向上。　画像解析・追尾処理部５３は、時系列の医用画像を画像処理して解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出する。力学モデル構築部５５は、時系列の形態指標と時系列の形状変形指標と時系列の医用画像とに基づいて解析対象領域の構造流体解析に関する力学モデルを暫定的に構築する。統計的同定部６１は、暫定的に構築された力学モデルに基づく血管形態指標の予測値及び血液流量指標の予測値の少なくとも一方が血管形態指標の観測値及び血液流量指標の観測値の少なくとも一方に整合するように潜在変数同定領域に関する潜在変数を同定する。算出部５７，５９は、潜在変数が潜在変数同定領域に割り当てられた力学モデルに構造解析、流体解析、あるいは、構造流体連成解析を施して時系列の力学的指標の予測値と時系列の血液流量指標の予測値との少なくとも一方を算出する。

Description

血管解析装置、医用画像診断装置、及び血管解析方法
　本実施形態は、血管解析装置、医用画像診断装置、及び血管解析方法に関する。
　三大疾患の一つである心疾患の原因となる冠動脈の狭窄や、脳動脈瘤、あるいはそれらの予兆となる頸動脈のプラークによる狭窄を、非侵襲あるいは低侵襲に予防・診断するための技術が望まれている。
　冠動脈の狭窄は、虚血性心疾患に至る重大な病変である。冠動脈の狭窄診断としては、カテーテルによる冠動脈造影検査（ＣＡＧ：Coronary Angiography）が主流である。冠動脈の器質的病変の診断指標として、心筋血流予備量比（ＦＦＲ：Fractional Flow Reserve）がある。ＦＦＲは、狭窄非存在下の最大冠血流に対する狭窄存在下の最大冠血流の比率として定義される。ＦＦＲは、狭窄近位部冠内圧に対する狭窄遠位部冠内圧の比率に略一致する。カテーテル先端に設けられた圧力センサが測定される。すなわち、ＦＦＲの測定は、カテーテル手術を必要とする。
　冠動脈の狭窄解析が心臓ＣＴで可能であれば、カテーテル手術によるＦＦＲの測定に比して、低侵襲、患者の負担低減、及び医療コストが節約できる。しかしながら、心臓ＣＴにおいては、ＣＴ画像に含まれるプラーク領域あるいは狭窄領域の大きさに基づく指標しか低侵襲に計測することしかできない。ＣＴ画像に基づいて狭窄前後の圧力差などを構造流体解析により計測できれば、狭窄（あるいはプラーク）が及ぼす影響の定量化が期待できる。
　冠循環の動態評価として臨床的には、超高速ＣＴ、シネアンギオグラム、超音波法、ＳＰＥＣＴ（シングルフォトンエミッショントモグラフィ）やＰＥＴ（ポジトロンエミッショントモグラフィ）を含む核医学イメージング、ＭＲＩ（核磁気共鳴画像法）などが開発導入され、診断や治療法の評価に役立っている。
　しかしながら、冠微小血管を医用画像診断装置で正確に捉えることは困難である。また、血管形状が鮮明であっても医用画像にノイズが含まれる場合や、生体組織の境界の閾値設定に曖昧性が存在する場合も多い。このように、医用画像診断装置から得られる血管形状は不確定性を有している。
　臨床応用で医用画像診断装置を活用する場合、冠微小血管より上流の大動脈起始部から冠動脈の太い領域のみを対象として解析が行われる場合も多い。冠動脈の血流が冠微小血管の緊張性（トーヌス）にも大きく影響をうけるため、太い領域の冠動脈の出口における流量あるいは圧力もしくはそれらの変化率といった流体解析の境界条件を適切に設定することが課題となる。また、冠動脈の血流は、心臓の拍動による機械的因子（拍動による全体的な動き、局所的な伸縮、ねじり、せん断変形による強制変位あるいは外力）を受ける。流体解析のみでは、心臓の拍動等の機械的因子の影響を考慮できないため、血流の流量分布や内圧分布を精確に計測することができない。一方、画像で捉えられる心臓および血管系を対象とし、機械的因子の影響を考慮した構造-流体連成解析も実施されている。しかし、構造－流体連成解析を行う場合でも、血液（造影剤を含む）の流体解析における血管の入口や出口の境界条件や血管やプラークの材料モデルを正しく設定することが困難な場合も多い。また、画像に描出されない微小血管が存在する場合には、微小血管が血流に与える影響を考慮できない場合もある。そのため、構造-流体連成解析の解析結果は、実際の血流や血管変形を再現できていない恐れがある。また、境界条件や負荷条件、材料モデルが適切でない場合や、血管が大きな動きを伴う場合、収束性や解析安定性に問題がある場合もある。このように、従来の血管の構造流体解析は、多大な解析リソースと解析時間とが必要となる場合や、解析結果の誤差が大きくなる場合があり、現実的に臨床の現場で活用するのに問題が生じる場合がある。
特開２００８－２４１４３２号公報
門岡ら（ITUジャーナル，心臓シミュレータが拓くテーラーメード医療～世界最先端の心臓シミュレータとその適用例のご紹介～，Vol.41,No.6）
　本実施形態に係る血管解析装置は、被検体の血管に関する時系列の医用画像のデータを記憶する記憶部と、前記時系列の医用画像に含まれる血管領域に解析対象領域を設定し、前記解析対象領域に潜在変数同定領域を設定する設定部と、前記時系列の医用画像を画像処理して前記解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出する算出部と、前記時系列の形態指標と前記時系列の形状変形指標と前記時系列の医用画像とに基づいて前記解析対象領域の構造流体解析に関する力学モデルを暫定的に構築する構築部と、前記暫定的に構築された力学モデルに基づく血管形態指標の予測値及び血液流量指標の予測値の少なくとも一方が予め計測された血管形態指標の観測値及び血液流量指標の観測値の少なくとも一方に整合するように前記潜在変数同定領域に関する潜在変数を同定する同定部と、を具備する。
　実施形態の目的は、血管（血液を含む）の構造流体解析および画像解析・追尾処理に関する精度の向上を可能とする血管解析装置、医用画像診断装置、及び血管解析方法を提供することにある。
図１は、本実施形態に係る医用画像診断装置（Ｘ線コンピュータ断層撮影装置）の概略的なブロック構成を示す図である。 図２は、本実施形態に係る構造流体解析の対象領域に関する力学モデルの一例を示す図である。 図３は、図１のシステム制御部の制御のもとに行われる構造流体解析処理の典型的な流れを示す図である。 図４は、図１の画像処理装置のブロック構成を示す図である。 図５は、血管の芯線に直交する断面を模式的に示す図である。 図６は、図４の画像解析・追尾処理による画像追尾処理に利用される血管芯線の形態の時間的変化を示す図である。 図７は、図４の画像解析・追尾処理による血管芯線の曲げ変形や回転変位の算出例を示す図である。 図８は、図４の画像解析・追尾処理による画像追尾処理を説明するための図であり、時刻ｔと時刻ｔ＋Δｔとの間における追尾処理の一例を示す図である。 図９は、図４の力学モデル構築部により構築される形状モデルの芯線に直交する断面を示す図である。 図１０は、図４の力学モデル構築部による形状モデルへの強制変位履歴の割り当てを説明するための図である。 図１１は、図４の力学モデル構築部による形状モデルへの強制変位履歴の他の割り当て方法を示す図であり、血管形状変形指標がねじれの場合における割り当て例を示す図である。 図１２は、図４の力学モデル構築部による形状モデルへの強制変位履歴の他の割り当て方法を示す図であり、血管形状変形指標が曲げの場合における割り当て例を示す図である。 図１３は、図４の統計的同定部による、階層ベイズモデル及びマルコフ連鎖モンテカルロ法による負荷条件（血管内の平均圧力）に関する事後分布算出と平均内圧の同定とを説明するための図である。 図１４は、図４の統計的同定部による、階層ベイズモデル及びマルコフ連鎖モンテカルロ法による材料モデルパラメータに関する事後分布算出と材料モデルパラメータ（血管壁の等価弾性率）の同定とを説明するための図である。 図１５は、図１の表示機器による力学的指標の一つである内圧の空間分布の表示例を示す図である。 図１６は、図１の表示機器による血液流量指標の一つである流速値の空間分布の表示例を示す図である。 図１７は、図１の表示機器により表示される左冠動脈起始部の血圧に関するグラフである。 図１８は、図１の表示機器により表示されるＬＣＸとＬＤＡとの分岐点付近の血圧に関するグラフである。 図１９は、図１の表示機器により表示される芯線方向に関する血圧変化に関するグラフである。 図２０は、図４の力学モデル構築部による強制変位履歴の他の割り当て例を示す図である。 図２１は、図４の力学モデル構築部による強制変位履歴の他の割り当て例を示す図である。 図２２は、図４の力学モデル構築部による強制変位履歴の他の割り当て例を示す図である。 図２３は、本実施形態に係る繊維グループの変形を示す図である。 図２４は、本実施形態に係る肉厚円筒の力学モデルの直交断面を示す図である。 図２５は、図２０の微小扇形要素の拡大図である。
以下、図面を参照しながら本実施形態に係わる血管解析装置、医用画像診断装置、及び血管解析方法を説明する。
　本実施形態に係る血管解析装置は、医用画像診断装置により発生された医用画像に含まれる血管領域を構造流体解析するためのコンピュータ装置である。本実施形態に係る血管解析装置は、医用画像診断装置に組み込まれていても良いし、医用画像診断装置とは別体のワークステーション等のコンピュータ装置であっても良い。以下、説明を具体的に行うため本実施形態に係る血管解析装置は、医用画像診断装置に組み込まれているものとする。
　本実施形態に係る医用画像診断装置は、被検体を撮像するための撮像機構を装備する如何なる種類の画像診断装置にも適用可能である。本実施形態に係る医用画像診断装置としては、例えば、Ｘ線コンピュータ断層撮影装置（Ｘ線ＣＴ装置）、磁気共鳴診断装置、超音波診断装置、ＳＰＥＣＴ装置、ＰＥＴ装置、放射線治療装置等が適宜利用可能である。以下、説明を具体的に行うため本実施形態に係る医用画像診断装置は、Ｘ線コンピュータ断層撮影装置であるものとする。
　図１は、本実施形態に係る医用画像診断装置（Ｘ線コンピュータ断層撮影装置）の概略的なブロック構成図である。図１に示すように、Ｘ線コンピュータ断層撮影装置は、ＣＴ架台１０とコンソール２０とを有する。ＣＴ架台１０は、コンソール２０の架台制御部２３からの制御に従ってＸ線で被検体の撮像部位を撮像する。撮像部位は、例えば、心臓である。ＣＴ架台１０は、Ｘ線管１１、Ｘ線検出器１３、及びデータ収集装置１５を有している。Ｘ線管１１とＸ線検出器１３とは、回転軸Ｚ回りに回転可能にＣＴ架台１０に装備されている。Ｘ線管１１は、造影剤が注入された被検体にＸ線を照射する。Ｘ線検出器１３は、Ｘ線管１１から発生され被検体を透過したＸ線を検出し、検出されたＸ線の強度に応じた電気信号を発生する。データ収集装置１５は、Ｘ線検出器１３から電気信号を読み出してデジタルデータに変換する。１ビュー毎のデジタルデータのセットは、生データセットと呼ばれている。複数のスキャン時刻に関する時系列の生データセットは、非接触データ伝送装置（図示しない）によりコンソール２０に伝送される。
　コンソール２０は、システム制御部２１を中枢として、架台制御部２３、再構成装置２５、画像処理装置２７、入力機器２９、表示機器３１、及び記憶装置３３を有している。
　架台制御部２３は、ユーザにより入力機器２９を介して設定されたスキャン条件に応じてコンソール２０内の各装置を制御する。
　再構成装置２５は、生データセットに基づいて被検体に関するＣＴ画像のデータを発生する。具体的には、まず、再構成装置２５は、生データセットに前処理を施して投影データセットを発生する。前処理としては、対数変換や不均一補正、キャリブレーション補正等が含まれる。次に、再構成装置２５は、投影データセットに画像再構成処理を施してＣＴ画像のデータを発生する。画像再構成アルゴリズムとしては、ＦＢＰ（filtered back projection）法等の解析学的画像再構成法や、ＭＬ－ＥＭ（maximum likelihood expectation maximization）法やＯＳ－ＥＭ（ordered subset expectation maximization）法等の逐次近似画像再構成等の既存のアルゴリズムが採用可能である。本実施形態において再構成装置２５は、時系列の投影データセットに基づいて時系列のＣＴ画像のデータを発生する。ＣＴ画像は、造影剤により造影された血管に関する画素領域（以下、血管領域と呼ぶことにする。）を含んでいる。なお、ＣＴ画像は、ＣＴ値の２次元空間分布を表現するスライスデータであっても良いし、ＣＴ値の３次元空間分布を表現するボリュームデータであっても良い。以下、ＣＴ画像はボリュームデータであるとする。時系列のＣＴ画像のデータは、記憶装置３３に記憶される。
　画像処理装置２７は、時系列のＣＴ画像に基づいて力学モデルを構築して構造流体解析を実行する。画像処理装置２７の処理の詳細については後述する。
　入力機器２９は、ユーザからの各種指令や情報入力を受け付ける。入力機器２９としては、キーボードやマウス、スイッチ等が利用可能である。
　表示機器３１は、ＣＴ画像や構造流体解析結果等の種々の情報を表示する。表示機器３１としては、例えばＣＲＴディスプレイや、液晶ディスプレイ、有機ＥＬディスプレイ、プラズマディスプレイ等が適宜利用可能である。
　記憶装置３３は、ハードディスク装置等の種々の記憶媒体により構成される。記憶装置３３は、時系列の投影データや時系列のＣＴ画像データ等の種々のデータを記憶する。例えば、記憶装置３３は、時系列のＣＴ画像データをＤＩＣＯＭ（digital imaging and communications in medicine）規格に準拠した医用画像ファイル形式で記憶する。また、記憶装置３３は、外部機器により収集された医用データを時系列のＣＴ画像データに医用画像ファイル内において関連付けて記憶しても良い。
　システム制御部２１は、中央演算処理装置（ＣＰＵ：central processing unit）や読み出し専用メモリ（ＲＯＭ：read only memory）、ランダムアクセスメモリ（ＲＡＭ：random access memory）を有する。システム制御部２１は、Ｘ線コンピュータ断層撮影装置の中枢として機能する。システム制御部２１は、ＲＯＭやＲＡＭに記憶されている血管解析プログラムを実行して本実施形態に係る血管構造解析処理を実行する。
　なお、システム制御部２１、画像処理装置２７、入力機器２９、表示機器３１、及び記憶装置３３は、血管解析装置５０を構成する。本実施形態のように血管解析装置５０は、医用画像診断装置（Ｘ線コンピュータ断層撮影装置）に組み込まれていても良いし、医用画像診断装置とは別体のコンピュータ装置であっても良い。血管解析装置５０が医用画像診断装置とは別体の場合、血管解析装置５０は、医用画像診断装置やＰＡＣＳ（picture archiving and communication systems）からネットワークを介して時系列のＣＴ画像等の医用データを収集すれば良い。
　次に、本実施形態に係る動作例を詳細に説明する。なお、本実施形態に係る血管解析装置、医用画像診断装置、及び血管解析方法は、心臓血管や頸動脈、脳動脈等の人体のあらゆる部位の血管を解析対象とすることができる。しかしながら、以下、説明を具体的に行うため本実施形態に係る解析対象は、心臓の血管であるとする。
　心臓の血管としては、例えば、冠動脈と大動脈とが挙げられる。冠動脈は、大動脈の冠動脈起始部から始まり心筋表面を走行し、心外膜側から内膜側に入り込む。冠動脈は、心筋の内膜において無数の毛細管に分岐する。分岐後、無数の毛細管は、再び統合して大心静脈を形成し、冠静脈洞に接続する。冠血管系は、他の臓器と異なり、心筋の収縮及び弛緩という力学的変化において灌流が保障されなければならないという特徴を有する。
　冠血流の特徴としては、心筋収縮による機械的血流阻害作用で冠動脈起始部の内圧が高くなる収縮期よりも、左心室拡張期に灌流圧が低下したときに多く流れることである。そのため、正常の冠動脈血流速波形は、収縮期と拡張期との二峰性であり、拡張期血流が優位である。肥大型心筋症や大動脈弁狭窄症では収縮期に逆行性波を認め、大動脈逆流症では収縮期順行波が大きくなるなど疾患によって特異的な血流波形を呈することが知られている。また、拡張期の順行性波形は左室拡張機能、特に左室弛緩と密接な関係がある。左室弛緩遅延例では拡張期波形のピークが後ろにずれ、また減速脚がゆるやかになる傾向がある。また、このような症例では、頻拍時には拡張期の冠血流は十分に増大できず、心筋虚血を助長すると考えられる。
　解剖学的に大動脈起始部から分岐する左右冠動脈に大動脈圧に等しい冠灌流圧、すなわち、冠動脈が分枝する大動脈起始部の圧力がかかることにより、冠血流が生じる。冠血流を決定するのは大動脈圧である駆動圧とともに冠血管抵抗が重要である。１４０～１８０μｍ以上の太い冠血管には冠血管低抗の２０％程度が存在するのに対し、１００～１５０μｍ以下の微小血管には抵抗成分の残りの多くが存在するといわれる。従って、いわゆる冠狭窄などのない場合には抵抗値は冠微小血管の緊張性（トーヌス）に左右される。
　血管抵抗因子としては、血管特性、動脈硬化、管狭窄、血液粘性、機械的因子があげられる。冠微小血管のトーヌスは、血管特性、心筋代謝（心筋酸素消費）、神経体液性因子、機械的因子、体液因子としての各種の血管作動性物質、血液粘性に規定され、さらに、心肥大、冠動脈硬化などを含めた様々な病変によっても影響され冠循環障害を起こす。
　冠動脈血流拍動は、冠動脈血流の拍動パターン、心筋収縮による心筋内血流の制御、機械的刺激に対する心筋内血管の反応に影響される。心筋収縮が血流を阻害する機序としては、心筋内圧の上昇、心筋内血管容量の変化、心筋内血管の圧迫が挙げられる。心筋拡張期の血流規定因子には、拡張期の冠動脈圧、拡張期の血管外力、心拍数、心周期に占める拡張期の割合、心筋弛緩が存在する。
　本実施形態に係る血管解析装置５０は、時系列のＣＴ画像に基づいて力学モデルを構築し、力学モデルを利用して心臓の血管についての構造流体解析を実行し、血管内の力学的指標や血管流量指標を高精度に算出する。精度の高い力学的指標や血管流量指標を算出するためには、力学モデルに精度の高い潜在変数を割り当てる必要がある。血管解析装置５０は、力学モデルを構築する際、初期的な力学モデルに逆解析を施して潜在変数を統計的に同定する。これにより血管解析装置５０は、高精度の潜在変数を決定することができる。力学的指標は、血管壁に関する力学的な指標を意味する。血管壁に関する力学的指標としては、例えば、血管壁の変位に関する指標、血管壁に生じる応力やひずみに関する指標、血管内腔に負荷される内圧分布に関する指標、血管の硬さなどを表す材料特性に関する指標等に分類される。血管の硬さなどを表す材料特性に関する指標は、血管組織の応力とひずみとの関係を示す曲線の平均的な傾き等が挙げられる。血液流量指標は、血管を流れる血液に関する血行動態の指標を意味する。血管流量指標としては、例えば、血液の流量、血液の流速、血液の粘性等が挙げられる。
　潜在変数は、例えば、血管の材料構成式、あるいは血液の材料構成式といった材料モデルのパラメータ（例えばヤング率やポアソン比など）、血管内腔に負荷される内圧分布などの負荷条件パラメータ、構造解析や流体解析の境界条件パラメータ、時系列の形態指標や形状変形指標の不確定性に関連するばらつき分布パラメータの少なくとも一つを含む。ここで、時系列の形態指標や形状変形指標の不確定性に関連するばらつき分布パラメータとは、医用画像データには、各ＣＴ値のノイズに起因したばらつき分布や、生体組織の境界閾値の曖昧性に起因する確率分布などが存在することを考慮して、種々の不確定性を確率分布として表現したものである。種々の不確定性としては、例えば、血管組織や血液の境界座標および特徴点（血管分岐部や造影剤分散配置など）の空間座標における不確定性、幾何学的構造パラメータ（芯線に垂直な断面の内腔半径など）の不確定性、医用画像データ自体（ＣＴ値や境界閾値など）の不確定性が挙げられる。
　力学モデルは、血管や血液の挙動を表現するための数値モデルである。力学モデルは、構造流体解析の手法に応じて異なるタイプを有している。例えば、力学モデルは、連続体力学モデルと簡易的力学モデルとに分類される。連続体力学モデルは、例えば、有限要素法（ＦＥＭ：finite element method）や境界要素法に用いられる。簡易的力学モデルは、例えば、材料力学に基づく材料力学モデルと流れ学に基づく流体力学モデルとに分類される。なお、以下の説明において特に言及しない場合、力学モデルのタイプについては特に限定しないものとする。初期的な力学モデルは、潜在変数の確率分布や変数範囲から得られる潜在変数のパラメータに関するサンプリング集合（各パラメータの組み合わせの集合）が割り当てられた力学モデルを意味するものとする。
　図２は、構造流体解析の対象領域（以下、解析対象領域と呼ぶ）に関する力学モデルＭ１の一例を示す図である。図２に示すように、力学モデルＭ１は、大動脈領域Ｒ１と右冠動脈領域Ｒ２と左冠動脈領域Ｒ３とを有している。血液は、大動脈から右冠動脈または左冠動脈へ流れる。
　図２に示すように、力学モデルＭ１の大動脈起始部側の末端が血流の入口に設定され、右冠動脈領域の末端と左冠動脈領域の末端とが血流の出口に設定される。入口と出口との各々に境界条件が設定される。入口に関する境界条件は、例えば、入口における血流の流速あるいは血流による圧力、もしくはそれらの変化率を含む。出口に関する境界条件は、例えば、出口における血流の流速あるいは血流による圧力、もしくはそれらの変化率を含む。大動脈、右冠動脈、及び左冠動脈の変形は、血流に起因する血管壁への力学的作用、心臓の拍動による血管壁への力学的作用（外力）、血管断面境界の負荷条件、血管壁の材料モデル、血管の無応力状態、及び血管壁の幾何学的形状等の様々な要因に依存する。ここで、血流に起因する血管壁への力学的作用は、例えば、血流に起因する内圧と血流に起因するせん断力とを含む。血流に起因する内圧により、血管半径方向あるは血管内腔面の垂直方向に変形する。心臓の拍動による血管壁への力学的作用と、血流に起因するせん断力は、血管芯線方向に関する伸縮やねじり、曲げ等の血管の変形は、負荷条件として、大動脈領域Ｒ１、右冠動脈領域Ｒ２、及び左冠動脈領域Ｒ３に割り当てられる。具体的には、血管芯線方向に関する伸縮やねじり、曲げ等の血管の変形は、強制変位（移動ベクトルや回転変位）あるいは荷重ベクトルの時間的変化により表現される。また、血流に起因する内圧に基づく血管半径方向あるいは内腔面の垂直方向に関する変形は、血管内腔への圧力分布の時間的変化として割り当てられる。
　大動脈領域Ｒ１、右冠動脈領域Ｒ２、及び左冠動脈領域Ｒ３には、構造流体解析において、強制変位による変位拘束条件が割り当てられる。これにより、構造流体解析における血管壁の変形自由度を縮小することができ、計算収束性の安定化や解析時間の短縮を実現できる。
　また、例えば、血管の形状の変形度合は、血管壁の材料に依存する。そのため、大動脈領域Ｒ１、右冠動脈領域Ｒ２、及び左冠動脈領域Ｒ３に材料モデルが割り当てられる。また、血管の形状の変形度合は、血管の無応力状態にも依存する。負荷条件の初期値として、血管の残留応力分布が割り当てられても良い。血管解析対象の数値計算用力学モデルの空間を離散化した節点集合および節点から構成される要素集合において、材料モデルや境界条件、負荷条件等の解析条件が同定される領域と、同定されない領域とに区分される。解析条件の同定領域内の節点には強制変位履歴の変位拘束条件が割り当てられ、材料モデルの同定領域には、血管壁表面（外表面）の節点のみに強制変位履歴の変位拘束条件が割り当てられ、血管壁内部には変位自由度を確保し変位拘束が割り当てられない。これにより、構造流体解析の変形自由度を抑制することができ、安定的かつ効率的に解析を行うことができる。ただし、血管壁表面上に緩衝用のダミーの要素集合を設け、その表面の節点に強制変位を与えてもよい。これにより、血管内腔にプラークなどによる突起がある場合や、血管分岐部など内圧による負荷が芯線方向の断面外の変形に影響を与える場合には、血管内腔だけではなく血管壁の形態指標も参照して血管に負荷される荷重ベクトルと内圧とを分離して同定することができる。また、生理学的には壁表面の脂肪層を模擬しており、一方、数値計算上は、血管壁表面に強制変位を与えることで血管壁内部に局所的に実際とは異なる高い応力が発生することを避ける効果もある。
　これら材料モデル、境界条件、及び負荷条件等の潜在変数に関するパラメータは、後述する力学モデルに基づく逆解析（統計的同定処理）により同定される。逆解析により同定された精確な潜在変数は、力学モデルに割り当てられる。精確な潜在変数が割り当てられた力学モデルにより、解析対象血管領域外の血管や心臓等の外部要因による当該解析対象血管領域への影響を加味した構造流体解析あるいは流体解析あるいは構造解析あるいは画像解析に基づく血行動態解析を実行することができる。血管解析装置５０は、力学モデルの構築に関し、逆解析による潜在変数の同定により、従来例に係る次の４点の困難を解決することができる。困難１．冠動脈の材料モデルの同定方法。困難２．心臓の形状の変形の冠動脈への影響の組み込み。困難３．冠動脈の境界条件の同定方法。困難４．医用画像データの不確定性に起因したばらつきを有する血管形状による画像解析や構造流体解析。この４点の困難の克服により、血管解析装置５０は、逆解析による潜在変数の同定を行わない従来の血管構造流体解析に比して、解析精度の向上を実現する。
　次に、本実施形態に係る構造流体解析処理の詳細について説明する。図３は、本実施形態に係るシステム制御部２１の制御のもとに行われる構造流体解析処理の典型的な流れを示す図である。図４は、画像処理装置２７のブロック構成を示す図である。
　図３に示すように、構造流体解析処理において、まず、システム制御部２１により記憶装置３３から処理対象の医用画像ファイルが読み出され、画像処理装置２７に供給される。医用画像ファイルは、時系列のＣＴ画像のデータの他に、当該被検体に関する血管内腔に関する圧力値のデータ、血液流量指標の観測値のデータ、及びプラーク指標を含んでいる。ＣＴ画像以外ではＭＲＩ画像や超音波エコー画像であってもよい。時系列のＣＴ画像のデータは、時系列のＣＴ値の３次元空間分布を表現するデータである。時系列のＣＴ画像は、例えば、１心拍で約２０枚、すなわち、約２０心位相分のＣＴ画像を含んでいる。
　図３に示すように、システム制御部２１は、画像処理装置２７に領域設定処理を行わせる（ステップＳ１）。ステップＳ１において画像処理装置２７の領域設定部５１は、時系列のＣＴ画像に含まれる血管領域に構造流体解析の解析対象領域を設定する。解析対象領域は、冠動脈に関する血管領域の任意の一部分に設定される。例えば、領域設定部５１は、ユーザによる入力機器２９を介した指示、または、画像処理により血管領域に解析対象領域と同定対象領域とを設定する。
　ここで、図５を参照しながら、血管の構造について説明する。図５は、血管の芯線に直交する断面（以下、血管断面と呼ぶ）を模式的に示す図である。図５に示すように、血管は、管状の血管壁を有している。血管壁の中心軸は芯線と呼ばれている。血管壁の内側は内腔と呼ばれている。内腔に血液が流れる。内腔と血管壁との境は血管内壁と呼ばれている。血管壁の外側には心筋等の血管周辺組織が分布している。血管壁と血管周辺組織との境は血管外壁と呼ばれている。血管壁内部にプラークが発生することがある。図５に示すように、プラークは、例えば、石灰化プラークや粥状プラーク等に分類される。粥状プラークは不安定プラークと呼ばれることもある。粥状プラークは、軟らかく、血管内壁が破裂して血栓として血管内部に染み出す危険性がある。従って、プラークの性状を把握することは臨床的に有用である。プラークの性状や存在領域は、医用画像ファイルに含まれるプラーク指標により特定可能である。プラーク指標は、例えば骨のＣＴ値を基準に正規化したＣＴ値の大きさなどにより相対的に判別することができる。しかしながら、血管内部のプラークの変形特性や硬さを解析するのは容易ではない。
　ステップＳ１が行われるとシステム制御部２１は、画像処理装置２７に画像解析・追尾処理を行わせる（ステップＳ２）。ステップＳ２において画像処理装置２７の画像解析・追尾処理部５３は、時系列のＣＴ画像に画像処理を施して時系列の血管形態指標と時系列の血管形状変形指標とを算出する。具体的には、画像解析・追尾処理部５３は、時系列のＣＴ画像に画像解析処理を施すことにより時系列の血管形態指標を算出し、時系列のＣＴ画像に追尾処理を施すことにより時系列の血管形状変形指標を算出する。
　より詳細には、画像解析・追尾処理部５３は、画像解析処理において、各ＣＴ画像から血管領域を抽出し、血管の内腔に関する画素領域（以下、血管内腔領域と呼ぶ）と血管壁に関する画素領域（以下、血管壁領域と呼ぶ）とを特定する。画像解析・追尾処理部５３は、血管形態指標として、血管の芯線に垂直な断面、あるいは血管内腔面に垂直な面が、血管内腔、血管壁、プラーク領域に交わる領域上の複数の画素の３次元座標を特定する。なお、血管形態指標は、３次元座標に限定されず、種々の幾何学的指標でも良い。例えば、芯線に垂直な断面における一定角度ごとの血管内腔の半径や直径及び０°の方向ベクトル、あるいは断面における全角度に対する平均面積や平均半径、あるいは、芯線方向に垂直な複数の断面で囲まれた血管内腔容積、あるいは内腔面に垂直な複数断面で囲まれた血管壁容積やプラーク容積等が本実施形態に係る幾何学的指標として利用可能である。
　追尾処理において、画像解析・追尾処理部５３は、ユーザからの入力機器２９を介した指示または画像処理により特徴点や特徴形状、代表点、画素等の複数の追跡点を設定する。追跡点は、具体的には、血管領域や血液領域、造影剤領域、プロトン領域等に設定される。例えば、画像解析・追尾処理部５３は、血管分岐部や表面の特徴形状などの追跡点集合を設定する。各時刻（各心位相）における画像解析・追尾処理部５３の追尾処理により得られた追跡点集合の変位データから、力学モデルの血管壁表面、血管壁内部、又は血管内腔における節点の変位の時間的変化を補間処理などにより算出し、算出された時間的変化を強制変位として割り当てる。また、例えば、画像解析・追尾処理部５３は、力学モデルに血管芯線上の節点を定義する。画像解析・追尾処理部５３は、力学モデルの血管壁表面、血管壁内部、又は血管内腔における節点の変位の時間的変化から、血管の芯線方向に関する伸縮やねじりや曲げに関する変形を抽出し、血管芯線と芯線とに垂直な断面における節点の強制変位として与えることで表現しても良い。このように、血管形状変形指標としては、力学モデルにおける各時刻の節点の強制変位データ（強制変位履歴）を特定する。
　以下、図６、図７、及び図８を参照しながら、画像解析・追尾処理を説明する。図６は、血管芯線の形態の時間的変化を示す図である。図６に示すように、例えば、時系列の医用画像は、１心拍につき２０枚のＣＴ画像を含んでいるものとする。すなわち、心位相０％から９５％まで５％おきにＣＴ画像が得られているものとする。画像解析・追尾処理部５３により各ＣＴ.画像から血管領域の芯線が抽出される。図６に示すように、芯線の形態は、心位相の経過に従って変化している。
　図８は、時刻ｔと時刻ｔ＋Δｔとの間における追尾処理の一例を示す図である。図８に示すように、血管芯線上にＰ１からＰ１０の力学モデルにおける節点が設定されており、各断面上の血管の力学モデルの節点と力学的につながっている。ただし、血液の力学モデルの節点とは独立である。血管の追跡点の変位データをもとに、血管芯線上のＰ１からＰ２０の節点の変位データを補間などの処理により算出し、各節点に強制変位が設定されているものとする。血管形状変形指標と血管形態指標について説明するため、節点Ｐ１３と節点Ｐ１４とにより規定される局所血管領域ＲＡを考える。時刻ｔにおいて、芯線方向に関する節点Ｐ１３と節点Ｐ１４との間の距離がＬであり、血管領域の半径がｒであるとする。画像解析・追尾処理部５３から、節点Ｐ１３と節点Ｐ１４の血管芯線方向の伸縮やねじりや曲げといった強制変位を抽出することにより、節点Ｐ１３の強制変位（３次元空間における移動変位と、芯線方向の回転変位）と節点Ｐ１４の強制変位（３次元空間における移動変位と、芯線方向の回転変位）を算出する。
　図７は、血管芯線の曲げ変形や回転変位の算出例を示す図である。図７に示すように、例えば、ねじり角は、線ａ及びｂから構成される面の法線方向ベクトルcの変化から算出しても良い。
　図７や図８に示すように、画像解析・追尾処理部５３は、追跡点の座標と移動ベクトルとに基づいて、芯線上の各節点の強制変位（３次元空間における移動変位と、芯線方向の回転変位）を算出し、血管形状変形指標を算出する。例えば、画像解析・追尾処理部５３は、隣り合う２つの節点の座標差の時間変化を芯線方向に関する伸縮距離ΔＬとして算出する。また、画像解析・追尾処理部５３は、各芯線上の節点について、当該節点と当該節点を含む血管領域断面上の他の節点（血管内腔あるいは血管壁あるいはプラーク領域における節点）との間の距離の時間変化を半径方向に関する伸縮距離Δｒとして算出する。また、画像解析・追尾処理部５３は、各追跡点について、当該追跡点の近傍の複数の追跡点の座標と移動ベクトルとに基づいて、芯線上の当該節点の芯線方向のねじれ角度Δθを算出する。
　また、画像解析・追尾処理部５３は、血液領域の造影剤やプロトンの画像追尾により、血液流量指標として、流速あるいは芯線方向断面の平均流速あるいは平均流量を算出してもよい。
　血管形状変形指標は、力学モデルにおける強制変位として利用される。以下、時系列の血管形態指標を形状履歴と呼び、時系列の血管形状変形指標を強制変位履歴と呼ぶことにする。
　ステップＳ２が行われるとシステム制御部２１は、画像処理装置２７に構築処理を行わせる（ステップＳ３）。ステップＳ３において画像処理装置２７の力学モデル構築部５５は、形状履歴（時系列の血管形態指標）と強制変位履歴（時系列の血管形状変形指標）と時系列の医用画像（ＣＴ画像，ＭＲＩ画像，超音波エコー画像などのＤＩＣＯＭデータ）とに基づいて解析対象領域に関する力学モデルを暫定的に構築する。力学モデルは、構造流体解析を行うための解析対象領域に関する数値モデルである。
　以下、ステップＳ３について詳細に説明する。力学モデル構築部５５は、まず、医用画像と形状履歴とに基づいて、力学モデル（数理モデル）を解くための形状モデルを構築する。形状モデルは、各時刻における血管領域の幾何学的構造を模式的に表現したものである。形状モデルは、例えば、複数の離散化領域に区分されている。各離散化領域の頂点は、節点と呼ばれる。力学モデル構築部５５は、時刻毎の医用画像に含まれる血管領域と血管形態指標とに基づいて時刻毎の形状モデルを構築しても良いし、特定の時相の医用画像に含まれる血管領域と血管形態指標とに基づいて時刻毎の形状モデルを構築しても良い。また、例えば、初期の負荷状態として、解析対象領域に対応する血管に残留応力が存在しないと仮定する場合、無応力状態の時相として、解析対象領域に対応する血管が最も収縮した時相を無応力状態であると仮定する。
　図９は、形状モデルの芯線に直交する断面を示す図である。図９に示すように、形状モデルは、芯線から外側に向けて血管内腔領域、血管壁領域を有している。プラークが存在する場合、血管壁領域にプラーク領域を設けても良い。また、血管周辺組織による血管への影響を考慮する場合、血管周辺組織のダミー要素を血管壁領域の外側に設けても良い。
　形状モデルが構築されると力学モデル構築部５５は、各潜在変数の確率分布や変数範囲から得られる潜在変数のパラメータに関するサンプリング値（例えばマルコフ連鎖モンテカルロ法などによる、各パラメータの組み合わせの集合からのサンプリング）を力学モデルに設定する。例えば、力学モデル構築部５５は、図２に示すように、大動脈領域Ｒ１の大動脈起始部側の末端に入口に関する境界条件の同定対象の領域（以下、境界条件同定領域）を設定し、右冠動脈領域Ｒ２の末端と左冠動脈領域Ｒ３の末端とに出口に関する境界条件同定領域を設定する。力学モデル構築部５５は、各境界条件同定領域に境界条件の確率分布や変数範囲から得られる境界条件のパラメータに関するサンプリング値を割り当てる。また、力学モデル構築部５５は、大動脈領域Ｒ１、右冠動脈領域Ｒ２、及び左冠動脈領域Ｒ３に材料モデルの同定対象の領域（以下、材料モデル同定領域と呼ぶ）及び負荷条件の同定対象の領域（以下、負荷条件同定領域と呼ぶ）を設定する。力学モデル構築部５５は、各材料モデル同定領域に材料モデルの確率分布や変数範囲から得られる材料モデルのパラメータに関するサンプリング値を割り当て、各負荷条件同定領域に負荷条件の確率分布や変数範囲から得られる負荷条件のパラメータに関するサンプリング値を割り当てる。血管は、流量が０でも残留応力が存在すると言われている。例えば、力学モデル構築部５５は、流量が０の場合の残留応力を負荷条件の初期値として解析対象領域に割り当てても良い。また、力学モデル構築部５５は、幾何学的構造に不確定性がある部位に、幾何学的構造の同定対象の領域（以下、幾何学的構造同定領域と呼ぶ）を設定しても良い。なお、幾何学的構造のパラメータは、幾何学的構造の不確定性に関連するばらつき分布パラメータ、あるいは医用画像データに内在するばらつき分布パラメータであり、各ＣＴ値のばらつき分布や生体組織の境界閾値のばらつき分布などであってよい。詳細は、後述するが、力学モデル構築部５５は、プラーク領域に材料モデルを設定しても良い。材料モデルの詳細については後述する。
　形状モデルが構築されると力学モデル構築部５５は、ステップＳ２において算出された時系列の血管形状変形指標、すなわち、強制変位履歴を形状モデルに割り当てる。潜在変数及び強制変位履歴が割り当てられた形状モデルを力学モデルと呼ぶことにする。
　図１０は、形状モデルＭ２は血管や血液の力学モデルの一部を示しており、力学モデルの節点への強制変位履歴の割り当てを説明するための図である。図１０は、形状モデルＭ２の一部分を示している。ただし、図１０では芯線がＭ２内に位置する場合を示しているが、芯線がＭ２外に位置する場合でも良い。図１０に示すように、形状モデルＭ２には複数の節点ＰＮ（ＰＮ１，ＰＮ２）が設定されている。芯線上の節点をＰＮ１と称し、血管や血液の力学モデルにおける節点をＰＮ２と称することにする。形状モデルＭ２は、ダミー要素表面、血管外壁、血管内壁、プラーク領域表面、プラーク領域内部、あるいは血液領域に設定される。力学モデル構築部５５は、形状モデルＭ２の各節点ＰＮ１に強制変位、すなわち、血管形状変化指標を時刻毎に割り当てる。
　具体的には、力学モデル構築部５５は、芯線上の隣り合う節点ＰＮ１と節点ＰＮ１とをビーム要素（あるいはリジッド要素）ＥＢで結ぶ。また、力学モデル構築部５５は、節点ＰＮ１と当該節点ＰＮ１を通る直交断面に含まれる他の節点ＰＮ２とをビーム要素ＥＢで結ぶ。力学モデル構築部５５は、節点ＰＮ１及びビーム要素ＥＢに各血管形状変形指標の形状変位方向に関する拘束条件を割り当てる。材料モデルや血管内腔の内圧を同定する領域では、強制変位としては、芯線方向に関する血管壁（あるいはダミー要素）表面の伸縮、血管壁（あるいはダミー要素）表面のねじれ、血管壁（あるいはダミー要素）表面の曲げ変形が挙げられる。例えば、材料モデルや血管内腔の内圧を同定しない領域では、芯線方向に関する強制変位だけではなく、半径方向に関する血管壁の時系列の伸縮（変位）も強制変位履歴として割り当てられる。また、血管内腔に突起がある場合や血管分岐部などのように、内圧が芯線方向断面外への変形に寄与する場合、当該領域には強制変位履歴が割り当てられずに、当該領域の周辺部（たとえばダミー要素の表面節点）のみに強制変位履歴が割り当てられる。また、力学モデル構築部５５、節点ＰＮ１及びビーム要素ＥＢに時系列の血管形状変形指標を強制変位履歴として割り当てられる。これにより、血管の全体及び局所に関する伸縮変形やねじれ変形や曲げ変形が表現される。
　なお、強制変位履歴の割当対象は、芯線上の節点及びビーム要素に限定されない。図１１及び図１２は、形状モデルへの強制変位履歴の他の割り当て方法を示す図である。図１１は、血管形状変形指標がねじれの場合における割り当て例を示し、図１２は、血管形状変形指標が曲げの場合における割り当て例を示している。図１１及び図１２に示すように、力学モデル構築部５５は、形状モデルの表面あるいは内部の節点に直接的に強制変位履歴を割り当てても良い。例えば、ステップＳ２において画像解析・追尾処理部５３は、特徴点についての伸縮量やねじり量や曲げ量等の血管形状変形指標を算出する。力学モデル構築部５５は、算出された血管形状変形指標を特徴点の周囲の節点へ補間（内挿や外挿）することにより直接的に割り当てる。
　本実施形態に係る画像処理装置２７は、ステップＳ３において暫定的に構築された力学モデルを用いて逆解析を施し、力学モデルに設定される潜在変数を統計的に同定する。統計的同定処理は、後述のステップＳ６において行われる。ステップＳ４及びＳ５は、それぞれ統計的同定処理に用いる血管形態指標及び血液流量指標を算出するために設けられている。
　ステップ３が行われるとシステム制御部２１は、画像処理装置に２７に血管応力解析処理を行わせる（ステップＳ４）。ステップＳ４において画像処理装置２７の血管応力解析部５７は、現段階の力学モデルに血管応力解析を施して時系列の血管形態指標の予測値を算出する。血管形態指標は、既述の血管形態指の何れであっても良いが、例えば、血管芯線方向に関する内腔領域の断面形状指標や血管壁の断面形状指標が用いられると良い。具体的には、内腔領域の断面形状指標は、内腔領域の注目画素の座標値、内腔領域の幾何学的構造パラメータ（内腔領域の半径、内腔領域の直径など）の少なくとも一つである。また、血管壁領域の断面形状指標は、具体的には、血管壁領域の注目画素の座標値、血管壁領域の幾何学的構造パラメータ（血管壁領域の半径、壁領域の直径など）の少なくとも一つである。なお、予測値は、力学モデルに血管応力解析を施して算出された血管形態指標の算出値を意味する。
　また、ステップ３が行われるとシステム制御部２１は、画像処理装置２７に血液流体解析処理を行わせる（ステップＳ５）。ステップＳ５において画像処理装置２７の血液流体解析部５９は、暫定的に構築された力学モデルに血液流体解析を施して時系列の血液流量指標の予測値を算出する。血液流量指標は、血流量または流速である。あるいは、血液流量指標は、血流量または流速の空間的または時間的な平均値であっても良い。なお、予測値は、力学モデルに血液流体解析を施して算出された血液流体指標の算出値を意味する。
　ステップＳ４及びＳ５が行われるとシステム制御部２１は、画像処理装置２７に同定処理を行わせる（ステップＳ６）。ステップＳ６において画像処理装置２７の統計的同定部６１は、ステップＳ４において算出された血管形態指標の予測値とステップＳ５において算出された血液流量指標の予測値とが、事前に収集された血管形態指標の観測値と血液流量指標の観測値とに整合するように力学的モデルの潜在変数のパラメータを統計的に同定する。
　図４に示すように、統計的同定部６１は、第１統計的同定部６１－１と第２統計的同定部６１－２とを有している。第１統計的同定部６１－１は、血管形態指標の予測値が血管形態指標の観測値に整合するように潜在変数のパラメータを統計的に同定する。第２統計的同定部６１－２は、血液流量指標の予測値が血液流量指標の観測値に整合するように潜在変数のパラメータを統計的に同定する。以下、第１統計的同定部６１－１と第２統計的同定部６１－２とを順番に説明する。
　具体的には、ステップＳ６において第１統計的同定部６１－１は、ステップＳ４において算出された血管形態指標の予測値と観測値とに基づくデータ分布を設定する。データ分布は、例えば、血管形態指標の予測値と観測値との誤差に関する多変量正規分布関数を示す。具体的には、統計的同定処理部６１－１は、力学モデルにおける各節点あるいは各要素について予測値と観測値との間の誤差に関する正規分布関数値を算出し、各正規分布関数値の積をデータ分布として設定する。データ分布は、時刻毎に個別に設定されても良いし、複数時刻まとめて設定されても良い。次に、第１統計的同定部６１－１は、力学モデルの潜在変数に事前分布（事前確率分布）を割り当てる。具体的には、材料モデル、境界条件、負荷条件、及び時系列の形態指標や形状変形指標の不確定性に関連するパラメータの各々に事前分布が割り当てられる。例えば、負荷条件のパラメータの一つである血管内腔に関する圧力値に関する事前分布が割り当てられる。圧力値の取り得る値の範囲（想定範囲）は、経験的に予め限定することができる。第１統計的同定部６１－１は、これら想定範囲内に限定して内圧値に関するモンテカルロシミュレーションを実行することにより、各離散化領域について内圧値の確率分布、すなわち、事前分布を算出する。また、力学モデル構築部５５は、事前分布として、芯線方向の圧力分布は滑らかであること、また、時間経過に伴う圧力変化が滑らかであること、血流の逆流がないことが観測されている場合には芯線方向の平均的な圧力変化の傾きが負であることを、例えば、多変量正規分布関数により数学的に表現した確率分布を事前分布として設定しても良い。これら想定範囲内に限定された確率分布に従って、負荷条件のパラメータに関するモンテカルロシミュレーションを実行することができ、力学モデルへ設定するための負荷条件（潜在変数）のサンプリング値を得ることができる。次に、第１統計的同定部６１－１は、各潜在変数について、事前分布とデータ分布とに統計的同定処理を施すことにより事後分布（事後確率分布）を算出する。統計的同定処理は、例えば、階層ベイズモデルやマルコフ連鎖モデルが挙げられる。そして、第１統計的同定部６１－１は、各潜在変数について、事後分布の最頻値や平均値等の統計値から各潜在変数のパラメータを同定する。例えば、上述の例の場合、血管内腔圧力値に関する事後分布が算出され、この事後分布から血管内腔圧力値の同定値が算出される。
　図１３は、階層ベイズモデル及びマルコフ連鎖モンテカルロ法による負荷条件（血管内の平均圧力）に関する事後分布算出と平均内圧の同定とを説明するための図である。図１３に示すように、血管起始部から延びる血管に石灰化プラーク領域と粥状プラーク領域とが設定されているものとする。石灰化プラーク領域は材料モデル同定領域に設定され、粥状プラーク領域に材料モデル同定領域に設定される。血管起始部から血管芯線方向に沿って進むにつれて血管内圧が降下する。血管芯線に沿って複数の節点が設定される。各節点を含む直交断面（節点断面）において内腔内圧の事後分布が算出され、事後分布の最頻値が同定される。
　なお、血管形態指標の観測値としては、例えば、ステップＳ２において算出された血管形態指標が用いられる。
　第２統計的同定部６１－２による処理は、データ分布の算出に用いる指標が異なるだけで第１統計的同定部６１－１による処理と同様である。すなわち、第２統計的同定部６１－２は、まず、ステップＳ５において算出された血液流量指標の予測値と観測値とに基づくデータ分布を設定する。次に、第２統計的同定部６１－２は、力学モデルの潜在変数に事前分布を割り当てる。例えば、血管に関する材料モデルのパラメータや血液に関する材料モデルのパラメータ、プラークに関する材料モデルのパラメータに関する事前分布が割り当てられる。これら材料モデルのパラメータとしては、例えば、弾性率などの材料モデルパラメータや、血液の構成式における粘性に関するパラメータが挙げられる。材料モデルのパラメータの想定範囲や確率分布は、経験的に予め設定することができる。第２統計的同定部６１－２は、各離散化領域について材料モデルのパラメータの確率分布、すなわち、事前分布を設定し、これら想定範囲内に限定した確率分布に従って、材料モデルのパラメータに関するモンテカルロシミュレーションを実行することができ、力学モデルへ設定するための材料モデルパラメータ（潜在変数）のサンプリング値を得ることができる。次に、第２統計的同定部６１－２は、各潜在変数について、事前分布とデータ分布とに統計的同定処理を施すことにより事後分布を算出し、算出された事後分布の統計値から各潜在変数のパラメータを同定する。例えば、上述の例の場合、材料モデルのパラメータに関する事後分布が算出され、この事後分布から材料モデルのパラメータの同定値が算出される。
　図１４は、階層ベイズモデル及びマルコフ連鎖モンテカルロ法による材料モデルパラメータに関する事後分布算出と材料モデルパラメータ（血管壁の等価弾性率）の同定とを説明するための図である。図１４に示すように、血管モデルは、図１３と同様であるとする。材料モデル同定領域に限定して、血管壁の材料モデルのパラメータ（例えば、等価弾性率）の事後分布が算出され、事後分布の最頻値が同定される。
　なお、血液流量指標の観測値は、例えば、大動脈に送り出される血流量変化であると仮定し、血管形態指標の観測値を、時系列のＣＴ画像から画像処理により計測される左心室の容積変化値（ＣＦＡ）として用いることができる。造影剤の冠動脈内注入後の造影剤の画像追尾により特徴点の移動量の時間的変化を算出することにより、流速や流量を算出してもよい。また、造影剤の血管芯線方向あるいは時間的な特定領域の濃度変化量を取得し、その濃度変化を各領域の芯線方向距離間隔で除した値や、濃度変化の時間的変化率から、流速や流量を算出してもよい。ＭＲＩの場合はプロトンの画像追尾を用い、超音波エコーの場合には、コントラストエコー図法などにより流量を算出する。
　また、解析対象領域の各画素の座標値が確定値であることを前提としない場合、すなわち、解析対象領域の幾何学的構造に不確定性があると仮定する場合、潜在変数に幾何学的構造を含めても良い。この場合、統計的同定部６１は、各節点の座標値の芯線方向に関する所定範囲内の変動、あるいは、解析対象領域の径の所定範囲内の変動を表現した正規分布などの確率分布を事前分布として設定すると良い。この場合、解析対象領域の形状が滑らかであること、また、芯線における節点の順番が不変であるという制約も事前分布として設定しても良い。これら想定範囲内に限定した確率分布に従って、幾何学的構造のパラメータに関するモンテカルロシミュレーションを実行することができ、力学モデルへ設定するための幾何学的構造の不確定性パラメータ（潜在変数）のサンプリング値を得ることができる。
　なお、各ステップＳ６において第１統計的同定部６１－１による統計的同定処理と第２統計的同定部６１－２による統計的同定処理との両方が行われなくても良い。すなわち、各ステップＳ６においては、第１統計的同定部６１－１による統計的同定処理と第２統計的同定部６１－２による統計的同定処理との何れか一方が行われても良い。
　また、上記の例においては、第１統計的同定部６１－１は、血管形態指標の予測値が血管形態指標の観測値に整合するように潜在変数のパラメータを統計的に同定し、第２統計的同定部６１－２は、血液流量指標の予測値が血液流量指標の観測値に整合するように潜在変数のパラメータを統計的に同定するとした。しかしながら、統計的同定部６１は、構造－流体連成解析に基づいて、血管形態指標の予測値と血液流量指標の予測値とが血管形態指標の観測値と血液流量指標の観測値とに整合するように潜在変数のパラメータを統計的に同定しても良い。統計的同定部６１による統計的同定処理のさらなる詳細については後述する。
　ステップＳ６が行われると画像処理装置２７の力学モデル構築部５５は、ステップＳ６において算出された潜在変数のパラメータを力学モデルに設定する。
　そしてシステム制御部２１は、同定終了条件が満たされたか否かを判定する（ステップＳ７）。ステップＳ７において同定終了条件が満たされていないと判定した場合（ステップＳ７：ＮＯ）、システム制御部２１は、ステップＳ４、Ｓ５、Ｓ６、及びＳ７を繰り返す。ここで、同定終了条件は、同定終了を判定するための指標（以下、同定終了指標と呼ぶ）が規定値に達するか否かにより表現される。同定終了指標としては、例えば、血管形態指標の予測値と観測値との差分値が挙げられる。この場合、システム制御部２１は、この差分値が既定値よりも大きい場合、同定終了条件が満たされていないと判定し、差分値が既定値よりも小さい場合、同定終了条件が満たされたと判定する。また、同定終了指標は、例えば、モンテカルロ法のサンプリング点の数でも良い。この場合、システム制御部２１は、このサンプリング点の数が既定値よりも小さい場合、同定終了条件が満たされていないと判定し、サンプリング点の数が既定値よりも大きい場合、同定終了条件が満たされたと判定する。同定終了条件が満たされた場合、力学モデル構築部５５は、その時点の最新の力学モデルを最終的な力学モデルに設定する。
　最終的な力学モデルが構築されると、力学モデル構築部５５は、血管形状変形指標の観測値、最終的な力学モデルにおける負荷条件のパラメータ、及び材料モデルのパラメータを関連付けたモデル（以下、関連モデルと呼ぶ）を算出する。関連モデルは、記憶装置３３に記憶される。関連モデルは、検索の容易性等のため患者情報や検査情報等に関連付けて記憶されると良い。なお、血管形態指標や血液流量指標の観測値、最終的な力学モデルにおける負荷条件のパラメータ、及び材料モデルのパラメータは、必ずしもモデルの形態で関連付けられる必要はなく、例えば、テーブルあるいはデータベースであっても良い。
　上記のステップＳ４、Ｓ５、Ｓ６、及びＳ７は、同一の同定法で反復しても良いし、異なる同定法で反復しても良い。異なる同定法で反復する場合、例えば、まず、簡易的力学モデルを利用して潜在変数を暫定的に同定し、次に、連続体力学モデルを利用して潜在変数を正確に同定しても良い。このように統計的同定処理を異なる手法で２段階に分けて行うことにより、潜在変数のパラメータを短時間で収束させることができる。簡易的力学モデルを利用する方法としては、内圧及び外圧を厚肉円筒の材料力学の式が挙げられる。また、簡易的力学モデルを利用する方法としては、ハーゲン・ポアズイユ流れ及び修正ベルヌーイの式を用いても良い。連続体力学モデルを利用する方法としては、ＦＥＭ構造流体解析が挙げられる。簡易的力学モデルを利用する同定法と連続体力学モデルを利用する同定法との詳細については後述する。
　ステップＳ７において同定終了受件が満たされたと判定した場合（ステップＳ７：ＹＥＳ）、システム制御部２１は、画像解析・追尾処理部５３に修正処理を行わせても良い（ステップＳ８）。ステップＳ８において画像解析・追尾処理部５３は、統計的同定法による逆解析で得られた潜在変数のもとで実施した構造流体解析結果（力学的指標の予測値及び血液流体指標の予測値）が観測値（力学的指標の観測値及び血液流体指標の観測値）に整合するように、時系列の医用画像に含まれる血管領域の形状を修正しても良い。表示機器３１は、修正後の時系列の医用画像に基づく診断結果を表示する。これにより、血管解析装置５０は、最終的な力学モデルを考慮した診断結果を表示することができる。あるいは、表示機器３１は、逆解析による同定とその構造流体解析により観測結果とが整合しない血管箇所・領域を画面に表示しても良い。例えば、血管の挙動の動きの速い心位相の画像はぼやけることが多く、医用画像に基づく画像解析により観測した血管形状には誤差が大きい箇所や領域が存在する。比較的、血管の挙動が安定した心位相における医用画像データはノイズが少ない。誤差分布の小さい血管形状データに基づいて、誤差が大きな心位相での血管形状を、力学モデルを用いることで、正しく内挿することができる。誤差が大きい血管箇所や領域について、正しく内挿した形状を、元データからのばらつき分布とともに、表示することができる。これにより、血管形状表示の安定性を確保できるとともに、形状の不確定性をユーザが認識できる。
　ステップＳ８が行われるとシステム制御部２１は、画像処理装置２７に血管応力解析処理を行わせる（ステップＳ９）。ステップＳ９において画像処理装置２７の血管応力解析部５７は、最終的な力学モデルに血管応力解析を施して時系列の力学的指標の予測値の空間分布を算出する。具体的には、離散化領域毎に力学的指標の予測値が算出される。
　また、ステップＳ８が行われるとシステム制御部２１は、画像処理装置２７に血液流体解析処理を行わせる（ステップＳ１０）。ステップＳ１０において画像処理装置２７の血液流体解析部５９は、暫定的に構築された力学モデルに血液流体解析を施して時系列の血液流量指標の予測値の空間分布を算出する。具体的には、離散化領域毎に血液流量指標の予測値が算出される。
　なお、力学的指標または血液流量指標として、ＦＦＲが算出されても良い。
　ステップＳ９及びＳ１０が行われるとシステム制御部２１は、表示機器３１に表示処理を行わせる（ステップＳ１１）。ステップＳ１１において表示機器３１は、ステップＳ９において算出された時系列の力学的指標の予測値とステップＳ１０において算出された時系列の血液流量指標の予測値とを表示する。例えば、表示機器３１は、時系列の力学的指標または時系列の血管流量指標を、時系列の力学モデルを当該予測値に応じた色で動画的に表示する。このため、表示機器３１は、各種の予測値とカラー値（例えば、ＲＧＢ）との関係を示すカラーテーブルを保持している。表示機器３１は、予測値に応じたカラー値をカラーテーブルを利用して特定し、特定されたカラー値に応じた色で当該予測値に対応する離散化領域を表示する。
　図１５は、力学的指標の一つである内圧の空間分布の表示例を示す図である。図１５に示すように、表示機器３１は、力学モデルを構成する各離散化領域を当該離散化領域に関する内圧値に応じた色で動画的に表示する。ユーザは、力学モデルを観察することにより、経時的且つ空間的に変化する力学的指標を色で把握することができる。
　図１６は、血液流量指標の一つである流速値の空間分布の表示例を示す図である。図１６に示すように、表示機器３１は、力学モデルを構成する各離散化領域を当該離散化領域に関する流速値に応じた色で動画的に表示する。ユーザは、力学モデルを観察することにより、経時的且つ空間的に変化する血液流量指標を色で把握することができる。
　例えば、血管内が完全に狭窄している場合、狭窄部位の内圧は、非狭窄部位の内圧よりも小さい。力学的指標として内圧が指定された場合、ユーザは、力学モデル上での局所的な色の違いにより狭窄の有無を判断することができる。また、狭窄部位の流量は、非狭窄部位の流用よりも小さい。血液流量指標として流量が指定された場合、ユーザは、力学モデル上での局所的な色の違いにより狭窄の有無を判断することができる。
　また、血管応力解析部５７は、プラーク領域の材料モデルパラメータの同定結果に基づいて力学的指標として硬さ値の空間分布を算出しても良い。この場合、表示機器３１は、プラーク領域に関する硬さ値の空間分布を力学モデル上において表示しても良い。また、表示機器３１は、プラーク領域周辺の内圧分布や応力分布、ひずみ分布を表示して良い。ユーザは、これらの表示をプラークの性状と破綻しやすさとを推定することに活用することができる。
　力学的指標及び血液流量指標の予測値は、力学モデルの離散化領域の色で表現する方法のみに限定されない。例えば、図１７、図１８、及び図１９に示すように、グラフで表示しても良い。なお、図１７は、左冠動脈起始部の血圧に関するグラフである。図１７のグラフの縦軸は正規化した血圧に規定され、横軸は心位相［％］に規定される。図１８は、ＬＣＸとＬＤＡとの分岐点付近の血圧に関するグラフである。図１８のグラフの縦軸は正規化した血圧に規定され、横軸は心位相［％］に規定される。図１９は、芯線方向に関する血圧変化に関するグラフである。図１９のグラフの縦軸は血圧比に規定され、横軸は大動脈からの距離［ｍｍ］に規定される。表示機器３１は、力学的指標及び血液流量指標の予測値をグラフで表示することにより、これら値をユーザに簡便に把握させることができる。
　ステップＳ１１が行われると構造流体解析処理が終了する。
　なお、図１０において、強制変位履歴は、形状モデルの芯線部と外壁部とに設定されるとしたが、強制変位履歴の設定箇所は、これに限定されない。例えば、強制変位履歴は、芯線部と外壁部との間の血管壁領域に設定されても良い。
　なお、強制変位履歴の拘束条件の割り当て対象は、境界条件及び材料モデルを同定するか否かに応じて切り分けられても良い。図２０、強制変位履歴の他の割り当て例を示す図であり、形状モデルの断面を示している。例えば、図２０の（ａ）に示すように、境界条件及び材料モデルを同定する場合、形状モデルの外壁部ＯＷ上の節点ＰＮ２のみに強制変位履歴が割り当てられ、血管壁領域ＲＶの節点ＰＮ３に強制変位履歴を割り当てないと良い。図２０の（ｂ）に示すように、境界条件及び材料モデルを同定しない場合、形状モデルの外壁部の節点ＰＮ２と血管壁領域ＲＶの節点ＰＮ３との両方に強制変位履歴が割り当てられると良い。この場合、芯線上の節点ＰＮ１に強制変位履歴が割り当てられる。また、外壁部ＯＷの節点ＰＮ２と節点ＰＮ１とをビーム要素ＥＢで結び、ビーム要素ＥＢ上の節点ＰＮ２及びＰＮ３にも強制変位履歴が割り当てられても良い。このとき、半径方向に関する収縮及び膨張は、ビーム要素ＥＢの伸縮変位で表現する。なお、血管内腔領域ＲＩには、強制変位履歴が割り当てられなくて良い。
　図２１は、強制変位履歴の割り当ての他の例を示す図であり、血管周辺組織のダミー要素ＲＤを含む形状モデルの断面を示している。図２１に示すように、ダミー要素ＲＤは、血管壁領域ＲＮの外側に設定される。形状モデルがダミー要素ＲＤを含む場合、血管壁領域ＲＮに加え、ダミー要素ＲＤにも節点ＰＮ４が設定される。節点ＰＮ４にも強制変位履歴が割り当てられる。力学モデル構築部５５は、境界条件及び材料モデルを同定する場合、血管壁領域ＲＶに含まれる節点ＰＮ３に強制変位履歴を割り当て、境界条件及び材料モデルを同定しない場合、強制変位履歴を割り当てなくても良い。節点ＰＮ３に強制変位履歴を割り当てる場合、内腔領域ＲＩに関する形状指標以外にも血管壁領域ＲＶに関する形状指標も参照して材料モデルの同定が行われる。
　図２２は、強制変位履歴の割り当ての他の例を示す図であり、プラーク領域ＲＰを含む形状モデルの断面を示している。図２２に示すように、プラーク領域ＲＰは、血管壁領域ＲＶに含まれる。プラーク領域ＲＰは、材料モデル同定領域に設定される。プラーク領域ＲＰについては、内腔形状指標、血管壁形状指標、及びプラーク指標を考慮して材料モデルが同定される。既述のように、プラーク指標は、例えば、超音波診断装置による組織性状診断により得られたプラークの性状に関するデータである。力学モデル構築部５５は、性状に応じてプラーク領域を複数の部分領域に区分し、複数の部分領域に個別に材料モデル同定領域を設定する。各部分領域には、予め、当該部分領域の性状に応じたパラメータ範囲が設定されると良い。既述の統計的同定処理により、各部分領域についての材料モデルパラメータが統計的同定部６１により同定される。そして、ステップＳ１０において表示機器３１が力学的指標として血管の材料特性に関する指標を表示することにより、ユーザは、プラークの性状を正確且つ容易に把握することができる。
　次に、潜在変数の一つである材料モデルについて詳細に説明する。血管の材料モデルとしては、弾性モデル、超弾性モデル、異方性超弾性モデル、粘性特性を考慮した超弾性モデルなどが適用可能である。異方性超弾性モデルとしては、例えば、Y.C.Funにより提案された数理モデルや、Holzapfel-Gasser構成式と呼ばれる数理モデルが適用できる。単位参照体積あたりのひずみエネルギーは以下の（１）式で表わされる。（１）式の第１項は、コラーゲンを含まない等方性基礎材料のせん断変形に関するエネルギーを表し、第２項は、コラーゲンを含まない等方性基礎材料の体積変形に関するエネルギーを表し、第３項は、コラーゲン繊維各グループの寄与（繊維方向の分散を考慮）を表している。
Figure JPOXMLDOC01-appb-M000001
　図２３は、繊維グループの変形を示す図である。図２３に示すように、円筒形上の外膜を仮定する。芯線方向ｚと周方向θとにより規定される面における平均方向Ａでの繊維グループの変形は、下記の（２）式で表される。
Figure JPOXMLDOC01-appb-M000002
　（１）及び（２）式中の材料モデルに関するパラメータとしては、以下の表１に示すように、材料パラメータと繊維分散パラメータとがある。材料パラメータとしては、Ｃ１０、Ｄ，Ｋ１、Ｋ２などが用いられ、繊維分散パラメータとしては、Ｋａｐｐａやγなどが用いられる。各パラメータのデフォルト値と制約条件とは表１に示す通りである。
Figure JPOXMLDOC01-appb-T000003
　血液の材料モデルは、以下の（３）式のようなＣａｓｓｏｎ構成式や（４）式のようなＨＢ構成式が好適である。
Figure JPOXMLDOC01-appb-M000004
Figure JPOXMLDOC01-appb-M000005
　これら材料モデルのパラメータは、上述のステップＳ６において統計的同定部６１により、血管形態指標及び血液流体指標を用いた統計的同定処理により同定される。
　次に、統計的同定部６１により行われる統計的同定処理の詳細について説明する。
　時系列の医用画像から計測される血管形態指標及び血液流量指標のような観測変数は、不確定性を有している。統計的同定部６１は、このような不確定性が存在する状況下における潜在変数の統計的同定法として、階層ベイズモデルとマルコフ連鎖モンテカルロ法とに基づく統計的手法を活用する。
　既述のように、ステップＳ６において統計的同定部６１は、ステップＳ４において算出された血管形態指標または血液流量指標の予測値と観測値とに基づくデータ分布を設定する。データ分布は、例えば、血管形態指標または血液流量指標の予測値と観測値との誤差に関する多変量正規分布関数を示す。データ分布は、時刻毎に個別に設定されても良いし、複数時刻まとめて設定されても良い。次に、統計的同定部６１は、形状モデルの強制変位と潜在変数とに事前分布を割り当てる。事前分布は、取り得る値の確率分布を示す。次に、統計的同定部６１は、潜在変数に関する数値シミュレーションのパラメータサーベイを実行し、潜在変数と血管形態指標または血液流量指標との関係を表現するモデルを構築する。例えば、材料モデルパラメータと内圧分布パラメータと血管形態指標または血液流量指標との関係がモデルに規定される。なお、血管形態指標または血液流量指標と潜在変数との関係は、モデルという形態ではなく、データベースあるいはテーブルにより規定されても良い。これらモデル、データベース、又はテーブルは、記憶装置３３に記憶される。このモデル、データベース、又はテーブルを利用して統計的同定部６１は、事前分布から血管形態指標または血液流量指標の確率分布を算出する。統計的同定部６１は、階層ベイズモデルとマルコフ連鎖モンテカルロ法とにより得られる事後分布から潜在変数を統計的に同定する。
　この同定問題は、具体的には、次の３条件を満足しない不良設定問題となる。３条件は、（１）解の存在が保証される、（２）解が唯一に定まる、（３）解がデータに対して連続的に変化し測定誤差に対して解が安定している、である。不良設定問題は、正規化理論とその拡張とからなる枠組み内で捉えると扱いやすい。不良設定問題を解くことは、標準正規化理論では不十分である。不良設定問題を解くためには、内部状態の不連続点を検出し、検出された不連続点を内部状態の推定に役立てる方法論が必要となる。この点で、本実施形態に係る統計的同定処理においてマルコフ確率場理論が有効となる。
　血管形態指標と血液流量指標とに不確定性が存在する環境下において統計的同定部６１は、適切な制約条件の下において潜在変数の確率分布パラメータを同定する。適切な制約条件を決定するためには、解の性質を事前に知っている必要がある。統計的同定部６１は、解空間の制約条件に関するデータベースをシミュレーションと観測値とに基づいて発生する。統計的同定部６１は、発生されたデータベースを利用して、超多自由度大規模問題に対してマルコフ確率場理論と階層ベイズモデルとに基づく統計的同定処理を実行する。制約条件となる事前分布の設定では、多くの数値実験結果に基づいてこれらの要因に関するパラメータの確率分布が並列的に個別に構成される。統計的同定部６１は、複数の確率分布を統合しデータの欠損を補間することで潜在変数のパラメータを同定する。この処理のために、統計的同定部６１は、マルコフ確率場理論を用いたモデルに基づく階層ベイズ法による推定を行う。解析対象とする構造の変形状態の実測結果から、同定した中間変数をもとに、ある負荷条件と境界条件における圧力や流量分布が推定できるという仕組みである。
　冠動脈の構造流体解析における材料モデル、境界条件、及び負荷条件の同定問題は、非線形逆解析と位置付けられ、解の一意性及び安定性が保証されない場合が多い。生体組織の材料特性や血圧の現実的に取り得る範囲は先験情報として想定できるため、これらを事前分布の確率分布として設定できる。また、圧力や変位は空間的時間的に滑らかであることも想定できるためこの情報も先験情報として事前分布の確率分布として設定できる。あるいは、血液の流れに逆流が生じていない事実を考慮できる場合は血管芯線方向の全体的な圧力分布の傾きは負（圧力降下が存在）であることも制約条件として用いて良い。負荷条件（内圧分布など）、境界条件、及び材料モデルに対して、時系列のＣＴ画像に基づく血管形状変形指標の観測値と力学モデルに基づく血管形状変形指標の予測値との２乗誤差分布をデータ分布として設定できる。観測可能な平均流量に関する２乗誤差分布もデータ分布として追加してもよい。これらの事前分布とデータ分布とに基づいて階層ベイズモデルとモンテカルロ法とを利用して事後分布を算出することができる。事後分布の発生確率や分散により、潜在変数のパラメータの同定値を得ることができる。発生確率が高く、分散が小さいほうが確信度合の高い同定値といえる。事後分布が多峰性分布となる場合でも、複数の同定値のうちの分散が小さい同定値を選択すれば良い。あるいは、複数の同定値が存在し得る場合、各同定条件で構造流体解析を実施し、それぞれの可能性を認識して、同定値や解析結果を診断や予防の指針情報として活用することができる。時系列のＣＴ画像にも誤差が含まれていることから、力学モデルの各節点に関する血管形態指標にも誤差が含まれる。このため、各血管形態指標を、例えば、時系列のＣＴ画像から計測された血管形態指標の予測値を平均値とした正規分布の確率変数として扱い、位置の空間的順序を保つという制約を含めた上で、事前分布を設定してもよい。また、潜在変数のパラメータの同定において、一意性がなく、複数の候補が考えられる場合がある。この場合、時系列のＣＴ画像から計測された血管形態指標の不確定性に従う乱数のサンプリング点に対する、潜在変数の同定値のサンプル集合の変動幅をチェックすることで、各同定値の候補のロバスト性（安定性）を判定する。各同定値の候補のロバスト性に基づいて最終的な同定値を決定しても良い。
　次に、力学モデルの詳細について説明する。既述のように、力学モデル構築部５５は、力学モデルの種類に応じて異なるタイプの力学モデルを構築することができる。連続体力学に基づくＦＥＭを用いる場合、力学モデル構築部５５は、血管壁の応力解析用のための形状モデル（ＦＥＭモデル）と血液の流体解析用のための形状モデル（ＦＥＭモデル）との両方を構築する。
　材料力学に基づく簡易的な同定法を用いる場合、材料力学における内圧を受ける厚肉円筒の式から近似的に圧力と弾性率と変位との関係を求める。この場合、芯線方向に配列された複数の離散化領域の各々について厚肉円筒近似が形状モデルとして用いられる。具体的には、力学モデル構築部５５は、芯線上に離散的に配列された節点を通る断面上の血管内腔形状と血管壁表面形状と断面中心とを特定する。次に力学モデル構築部５５は、血管内腔形状と血管壁表面形状とに基づいて平均面積、内腔の平均半径、及び平均壁厚を算出する。そして力学モデル構築部５５は、平均面積、内腔の平均半径、及び平均壁厚に基づいて、各離散化領域の血管領域に厚肉円筒近似を施して形状モデルを構築する。
　流れ力学に基づく簡易的な同定法を用いる場合、血流の平均圧力と平均流量とを近似的に求めるため、流体力学における修正ベルヌーイの式、あるいは、Hagen-Poiseuille Flow（ハーゲン・ポアズイユ流れ）の式を用いる。この場合、複数の離散化領域の各々の血圧差と流量との関係を近似的に求めるための形状モデルが構築される。ここで、血圧差は、入口の圧力と出口の圧力との圧力差であり、流量は、単位時間あたり入口流量（あるいは流速）と出口流量（あるいは流速）とを意味する。ただし、力学モデル構築部５５は、ステップＳ２において算出された芯線方向に関する伸縮距離とねじれ角とに基づいて、各節点に関する移動ベクトルと当該節点に隣接する節点に関する回転変位とを各時刻で対応させると良い。
　次に、連続体力学モデルを用いる構造流体解析について説明する。
　血管及び血液（血管と血液を構成する物質）についての運動学は、その運動を引き起こしている力とは無関係である。血管及び血液についての運動学の基本概念は、位置、時間、物体、運動、及び変形し得る物質の集合についての直感的概念を数学的用語へ抽象化したものである。血管及び血液に関する変形及び運動の局所的解析を支配する基礎的運動学テンソルは、変形勾配テンソルFと速度勾配テンソルLとである。変形勾配テンソルFは、運動する血管及び血液の物質要素に生じる大きさ及び形の変化を規定する。変形勾配テンソルFは、回転テンソル（正格直交テンソル）Rとストレッチングテンソル（正値対称テンソル）U,Vとの積で表わされる。ストレッチングテンソルU,Vは、まず、基本形態での正規直交ベクトルRによって定められる方向に関するストレッチを課し、次に、正規直交ベクトルRによって与えられる剛体回転を課すことによりもたらされる。なお、ストレッチを課す順番と剛体回転を課す順番とは逆でも良い。速度勾配テンソルLは、基準形態に依存せず現在形態にのみ依存する。速度勾配テンソルLは、運動する血管及び血液の物質要素に生じる大きさ及び形の変化が生じる速度を規定する。速度勾配テンソルLは、ひずみ速度テンソルD（対称テンソル）とスピンテンソル（反対称テンソル）とに分離できる。ひずみ速度テンソルDは、物体がその現在形態をちょうど通過するときのストレッチの変化率を表す。スピンテンソルは、物体がその現在形態をちょうど通過するときの回転の変化率を表す。血管の力学モデルの外表面節点（あるいは積分点）の一部に関する変形勾配テンソルと速度勾配テンソルとに変位拘束条件（時間的変化を含む）を割り当てる。力学モデルの内腔節点（あるいは積分点）に関する予測値（変形勾配テンソル、速度勾配テンソル、あるいはそれらの関数値（例えば、変位や面積でも良い））と、観測値（観測データから得られた変形勾配テンソル、速度勾配テンソル、あるいはそれらの関数値）とが整合するように、材料モデルのパラメータ、負荷条件（力学モデルの内腔における表面力ベクトル）、及び境界条件（血管境界における力ベクトル）を同定する。ここで、血管内部の応力の初期状態は、予め仮定しておいてもよいし、同定してもよい。
　連続体力学に基づく力学モデルは、運動中の血管および血液における質量、運動学、運動量、角運動量、及びエネルギーの平衡を表す方程式を基礎としている。質量、力、熱、及び内部エネルギーという概念が基本である。平衡則とは、運動量の時間全微分が物体力と接触力との和に等しくなり、また、角運動量の時間全微分が物体トルクと接触トルクとの和に等しくなり、また、運動エネルギーと内部エネルギーとの時間的変化が、仕事率（力学的エネルギー）と単位時間あたりの熱供給と熱流束との和に等しくなる、ということを意味している。平衡則、構成式、及び跳躍条件から、変形勾配テンソル、速度勾配テンソル、及び応力テンソルといった場の方程式を導き、血管及び血液の力学モデルを記述できる。また、血管及び血液におけるひずみ場は、適合条件を満たす。ここで、構成式とは、密度、内部エネルギー、速度ベクトル、応力テンソル、熱流束ベクトル、及び温度からなる１０個のスカラ方程式の組の関係を与えるものである。場の釣り合い式の中の１７個のスカラ場、つまり、密度、内部エネルギー、速度ベクトル、応力テンソル、及び熱流束ベクトルのうち、場の釣り合いにより８個のスカラ関係が与えられ、温度を加えた残りの未知量の関係を与える。ただし、物体力ｂと熱源ｒとは既知としている。これらのスカラ場を与える方程式のパラメータが材料モデルパラメータである。
　連続体力学に基づく力学モデルは、有限要素法又は境界要素法による数値解析法により変位ベクトル、応力テンソル、ひずみテンソル、及び速度ベクトル等の場の近似解を、与えられた境界条件、負荷条件、及び材料モデルのもとに算出することができる。
　構造－流体連成解析において、構造及び流体の方程式を解く方法は、一体型解法（monolithic method）と分離型解法（partitioned method）とのいずれでも良い。また、構造と流体との間の境界面での連成は、弱連成でも強連成でもよい。また、流体解析においては、ＡＬＥ法に代表されるような境界面追跡型の手法で移動境界を扱ってもよいし、Immersed Boundary method、Immersed Finite Element Method、Fictitious Domain Methodなどの境界面補足型の手法でもよい。
　次に、簡易的な力学モデルの例として、内圧と外圧とを受ける厚肉円筒の材料力学の式と、ハーゲン・ポアズイユ流れ及び修正ベルヌーイの式とについて詳細に説明する。
　まず、図２４と図２５とを参照しながら、厚肉円筒の材料力学の式について説明する。図２４は、肉厚円筒の力学モデルの直交断面を示す図であり、図２５は、図２４の微小扇形要素の拡大図である。内半径ra、外半径rbの厚肉円筒に内圧paと外圧pbとが作用する場合の応力やひずみ、変位などの式を説明する。E及びνは材料モデルパラメータである。Eは弾性率、νはポアソン比を表している。厚肉円筒では，半径応力σrも考慮し、円周応力σθの半径方向分布も考慮する必要がある。以降では、軸方向のひずみεzは、断面の位置及び向きに関して一様であるとする。また円筒断面は、軸対称である場合について説明するが任意形状であっても良い。円筒断面が軸対称である場合、平衡条件は、任意の断面の半径方向についてのみ考えれば良い。任意の断面上で半径rとr+drの同心円筒と中心角dθで切り取られた単位厚さ1の微小扇形要素について、半径方向に関する力の平衡を考える。変形も軸対称であるため、ab面及びbc面には、せん断応力が生じないので、垂直応力のみ作用する。従って、半径方向に関する力の平衡は、以下の（５）式のように表現することができる。
σr rdθ+2σθr sin(dθ/2)-（σr+(dσr/dr)dr）(r+dr)dθ=0　　…（５）
　ここで、drはrより小さく、dσrはσrより小さいので、（５）式に含まれる高次の微小項を省略し、sin(dθ/2)≒dθ/2とすれば、（５）式は、以下の（６）式のように表現することができる。
rdσr/dr+σr-σθ=0　　…（６）
　半径rにおける半径方向の変位をuとすれば、u+drでの同方向の変位はu+(du/dr)drとなるので、半径rにおける半径方向のひずみεrはεr=du/drとなる。また、半径方向の変位uによって、半径rの円は半径r+uの円になる。従って、円周ひずみεθは、以下の（７）式のように表現することができる。
εθ=(2π(r+u)-2πr)/2πr=u/r　　…（７）
　また、応力とひずみの関係式から、以下の（８）式または（９）式が得られる。
d2u/dr2+(1/r)(du/dr)-u/r2=0　　…（８）　
d2u/dr2+d(u/r)/dr=0　　…（９）
　（８）式または（９）式を積分すると以下の（１０）式を得ることができる。
u=c1r+c2/r　　…（１０）
　これにより以下の（１１）式、（１２）式、（１３）式が得られる。
　σr=(E/((1+ν)(1-2ν)))(c1-(1-2ν)(c2/r2)+νεz) 　　…（１１）　
　σθ=(E/((1+ν)(1-2ν)))(c1+(1-2ν)(c2/r2)+νεz) 　　…（１２）　
　σz=(Eν/((1+ν)(1-2ν)))(2c1+((1-ν)/ν)εz) 　　…（１３）
　（１１）式、（１２）式、（１３）式における定数c1，c2は、周辺条件、すなわち、円筒の内周r=raでσr=-pa、外周r=rbでσr=-pbから定めることができる。この周辺条件により、（１１）式、（１２）式、及び（１３）式から以下の（１４）式、（１５）式、（１６）式をそれぞれ得ることができる。また、変位uは、以下の（１７）式のように表現することができる。
　σr=(1/(rb 2-ra 2))(ra 2(1-rb 2/r2)pa-rb 2(1-ra 2/r2)pb) 　　…（１４）　
　σθ=(1/(rb 2-ra 2))(ra 2(1+rb 2/r2)pa-rb 2(1+ra 2/r2)pb) 　　…（１５）　
　σz=2ν(ra 2pa-rb 2pb)/(rb 2-ra 2)+Eεz=ν(σr+σθ)+Eεz　　…（１６）　
　u=((1+ν)(1-2ν)/E)((ra 2pa-rb 2pb)/(rb 2-ra 2))r　
　　　 +((1+ν)/E)((ra 2rb 2)/((rb 2-ra 2)r))(pa-pb)-νεzr　　…（１７）
　上記の（１６）式及び（１７）式は、εzの項が含んでいる。従って、（１６）式は、対象とする円筒の境界における拘束条件に応じて異なる。例えば、対象とする円筒の両端が拘束される場合、εz=0となり、両端が開放され、σz=0となる。
　rθ面内のせん断応力τrは、以下の（１８）式のように表現することができる。
　τr=(1/2)(σr-σθ)=((ra 2rb 2)/((rb 2-ra 2)r2))pb　　…（１８）
　θz面内のせん断応力τr’は、σz=0のとき最大になるのでσz=0とすると、以下の（１９）式のように表現することができる。
τr’=(1/2)|σθ|=(1/2)((rb 2)/(rb 2-ra 2))(1+ra 2/r2)pb　　…（１９）
　次に、ハーゲン・ポアズイユ流れ及び修正ベルヌーイの式について説明する。
　流体が円管内へ流入すると、下流に進むにつれて圧力は降下するとともに、流れの速度分布も徐々に変化する。流れが管内へ流入すると管壁から境界層が発達し、下流へ進むにつれて境界層は厚さを増し、ついには管内の流れは境界層に覆われる。このため、速度分布は管入口のほぼ平らな分布から下流の放物形分布へと変化し、それ以降の速度分布は変化しなくなる。この状態を完全に発達した流れと呼び、管摩擦損失による圧力降下も一定の割合になる。流れが管入口から発達した流れに達する区間を助走区間または入口区間といい、その区間の長さを助走距離または入口長さと呼ぶ。完全に発達した流れの速度分布は下流方向へ変化しないから、以下の（２０）式に示すように、管摩擦損失によって生じる圧力損失ΔPの作用力と流体の粘性によって生じるせん断応力τの摩擦力は釣り合うことになる。円管内流れの場合、レイノルズ数Reがおよそ2300以下のとき層流となるといわれている。
Figure JPOXMLDOC01-appb-M000006
　速度分布uは層流の場合、以下の（２１）式のように軸対称な回転放物面で表わされる。
Figure JPOXMLDOC01-appb-M000007
　流量Qは、速度分布uを管断面全体にわたって積分することにより、以下の（２２）式のように表現することができる。
Figure JPOXMLDOC01-appb-M000008
断面平均流速vは、以下の（２３）式のように表現することができる。
Figure JPOXMLDOC01-appb-M000009
　圧力勾配は、流れ方向へ一定で、圧力は減少する。圧力勾配は、管長lの圧力降下Δpを用いると、以下の（２４）式のように表現できる。（２４）式に（２２）式を代入すると（２５）式が得られる。（２５）式は、流量Qが圧力損失Δpに比例することを意味する。この関係を満たす流れをハーゲン・ポアズイユ流れという。
Figure JPOXMLDOC01-appb-M000010
　損失がある場合の修正ベルヌーイの式は、以下の（２６）式のように表現することができる。
Figure JPOXMLDOC01-appb-M000011
　Δｈは損失ヘッドを表している。損失ヘッドは、血管内壁の摩擦損失ヘッド、狭窄部（広がり管や細まり管の損失ヘッド）や分岐部、曲り管部などの損失ヘッド、血管の解析対象領域入口／出口部の入口／出口損失ヘッド、管内非定常流れの圧力変動に起因した水撃による損失ヘッドのいずれかがある。例えば、管摩擦損失ヘッドは以下の式で表わされる。
Figure JPOXMLDOC01-appb-M000012
　ここで、lは管の長さ、dは管内径、vは管内平均速度、λは管摩擦係数である。管摩擦係数λは、流れが層流の場合にはレイノルズ数Reによって、乱流の場合にはレイノルズ数Reと表面粗さとによって定まる値になる。流体の粘性によって摩擦抵抗が必ず作用する。この摩擦抵抗は、流れを駆動する動力あるいはエネルギーを消費することになるのでエネルギー損失になる。
　ハーゲンポワズイユ流れの式をダルシー・ワイズバッハの式のように変形すると、圧力差と流速および管内径とに関する以下の（２８）式が得られる。
Figure JPOXMLDOC01-appb-M000013
　以上でハーゲン・ポアズイユ流れ及び修正ベルヌーイの式の説明を終了する。
　上述のように、画像処理装置２７は、この材料力学の式とハーゲン・ポアズイユ流れ及び修正ベルヌーイの式とを利用して潜在変数を同定することができる。例えば、力学モデルとして、血管の変形を厚肉管の材料力学の式を用い、管径変化を内圧変化と弾性率とにより表現する場合について考える。無応力状態を初期形状（例えば、血管が最も収縮する状態）と仮定した場合、血管壁及びプラークの弾性率をある値に設定すると、血管内腔の平均半径などの血管形状変形指標の観測値の時間的変化量と内圧の変化量との関係式が得られる。血管形状変形指標の観測値は、時系列のＣＴ画像から計測される。この血管形状変形指標の観測値の時間的変化量に合致するように血管の内圧分布の時間的変化が決定される。この内圧分布の下に血液の流体解析を行うことで血管流量指標の予測値が計測される。この血管流量指標の予測値が観測値に一致しない場合、画像処理装置２７は、最初に決めた血管壁あるいはプラークの弾性率を変更して、さらに同様の解析を行う。これを繰り返すことにより、画像処理装置２７は、血管形状変形指標の観測値と血液流量指標の観測値とに整合する血管壁およびプラークの弾性率、内圧分布、流体解析の圧力境界条件などの潜在変数を決定することができる。この決定方法をより効率的かつ安定的に行うために、階層ベイズモデルとマルコフ連鎖モンテカルロ法とによる統計的同定手法を用いてもよい。
　上記のように、本実施形態に係る血管解析装置５０は、記憶装置３３、領域設定部５１、画像解析・追尾処理部５３、力学モデル構築部５５、及び統計的同定部６１を有している。記憶装置３３は、被検体の血管に関する時系列の医用画像のデータを記憶する。領域設定部５１は、時系列の医用画像に含まれる血管領域に解析対象領域を設定する。画像解析・追尾処理部５３は、時系列の医用画像を画像処理して解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出する。力学モデル構築部５５は、時系列の形態指標と時系列の形状変形指標と時系列の医用画像とに基づいて解析対象領域の構造流体解析に関する力学モデルを暫定的に構築する。統計的同定部６１は、暫定的に構築された力学モデルに基づく血管形態指標及び血液流量指標の予測値が予め計測された血管形態指標及び血液流量指標の観測値に整合するように解析対象領域に関する力学モデルの潜在変数を同定する。
　上記構成により、本実施形態に係る血管解析装置５０は、材料モデル、境界条件、負荷条件、及び幾何学的構造等の潜在変数を血管形状変形指標と血液流量指標とを用いた逆解析により同定することができる。血管解析装置５０は、潜在変数を変更しながら逆解析を反復して行うことにより、上述の４点の困難、すなわち、１．冠動脈の材料モデルの同定方法、２．心臓の形状の変形の冠動脈への影響の組み込み、３．冠動脈の境界条件の同定方法、４．不確定性を有する血管形状を利用した材料モデルや負荷条件、境界条件の同定方法、を全て加味した潜在変数を同定することができる。従って、血管解析装置５０は、ＣＴ画像に描出されない血管や心臓等の外部要因による影響を加味した構造流体解析を実行することができる。
　かくして、本実施形態によれば、血管の構造流体解析の精度の向上が可能となる。
　本発明のいくつかの実施形態を説明したが、これらの実施形態は、例として提示したものであり、発明の範囲を限定することは意図していない。これら新規な実施形態は、その他の様々な形態で実施されることが可能であり、発明の要旨を逸脱しない範囲で、種々の省略、置き換え、変更を行うことができる。これら実施形態やその変形は、発明の範囲や要旨に含まれるとともに、特許請求の範囲に記載された発明とその均等の範囲に含まれるものである。
　１０…ＣＴ架台、１１…Ｘ線管、１３…Ｘ線検出器、１５…データ収集装置、２０…コンソール、２１…システム制御部、２３…架台制御部、２５…再構成装置、２７…画像処理装置、２９…入力機器、３１…表示機器、３３…記憶装置、５０…血管解析装置、５１…領域設定部、５３…画像解析・追尾処理部、５５…力学モデル構築部、５７…血管応力解析部、５９…血液流体解析部、６１…統計的同定部

Claims (15)

  1. 　被検体の血管に関する時系列の医用画像のデータを記憶する記憶部と、
    　前記時系列の医用画像に含まれる血管領域に解析対象領域を設定し、前記解析対象領域に潜在変数同定領域を設定する設定部と、
    　前記時系列の医用画像を画像処理して前記解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出する算出部と、
    　前記時系列の形態指標と前記時系列の形状変形指標と前記時系列の医用画像とに基づいて前記解析対象領域の構造流体解析に関する力学モデルを暫定的に構築する構築部と、
    　前記暫定的に構築された力学モデルに基づく血管形態指標の予測値及び血液流量指標の予測値の少なくとも一方が予め計測された血管形態指標の観測値及び血液流量指標の観測値の少なくとも一方に整合するように前記潜在変数同定領域に関する潜在変数を同定する同定部と、
    　を具備する血管解析装置。
  2. 　前記同定部は、前記潜在変数として、前記潜在変数同定領域に関する材料モデル、前記潜在変数同定領域の血液流入口及び血液流出口に関する境界条件、前記潜在変数同定領域の血流による内圧分布を含む負荷条件、及び前記潜在変数同定領域の幾何学的構造の少なくとも一つを同定する、請求項１記載の血管解析装置。
  3. 　前記同定部は、
    　前記暫定的に構築された力学モデルに応力解析を施して前記血管形態指標の予測値を推定する血管応力解析部と、
    　前記暫定的に構築された力学モデルに流体解析を施して前記血液流量指標の予測値を推定する血液流体解析部と、
    　前記血管形態指標の予測値が前記血管形態指標の観測値に整合するように前記潜在変数を統計的に同定する第１統計的同定部と、
    　前記血液流量指標の予測値が前記液流量指標の観測値に整合するように前記潜在変数を統計的に同定する第２統計的同定部と、
    　を有する、請求項２記載の血管解析装置。
  4. 　前記構築部と前記同定部とを制御する制御部をさらに備え、
    　前記制御部は、前記同定部により潜在変数が同定される毎に所定の同定終了指標が既定値に到達しているか否かを判定し、前記同定終了指標が前記既定値に到達していないと判定した場合、前記構築部を制御して前記同定された潜在変数を利用して力学モデルを再構築させ、前記同定部を制御して前記再構築された力学モデルの潜在変数を同定させ、前記同定終了指標が前記既定値に到達した判定した場合、前記構築部を制御して最新の力学モデルを最終的な力学モデルに設定させる、
    　請求項１記載の血管解析装置。
  5. 　前記最終的な力学モデルに応力解析を施して時系列の力学的指標を算出する血管応力解析部と、
    　前記時系列の力学的指標を表示する表示部と、
    　をさらに備える請求項４記載の血管解析装置。
  6. 　前記最終的な力学モデルに流体解析を施して時系列の血液流体指標を算出する血液流体解析部と、
    　前記時系列の血液流体指標を表示する表示部と、
    　をさらに備える請求項４記載の血管解析装置。
  7. 　前記最終的な力学モデルに基づく力学的指標と血液流量指標との少なくとも一方の予測値が力学的指標と血液流量指標との少なくとも一方の観測値に整合するように前記時系列の医用画像に含まれる血管領域を修正する修正部、をさらに備える請求項４記載の血管解析装置。
  8. 　前記構築部は、前記時系列の形態指標と前記時系列の医用画像とに基づいて前記解析対象領域の形状モデルを構築し、前記形状モデルに前記時系列の形状変形指標を強制変位として割り当て、前記潜在変数を割り当てることにより前記力学モデルを構築する、請求項２記載の血管解析装置。
  9. 　前記構築部は、前記力学モデルを材料モデル同定領域と非同定領域とに区分し、前記材料モデル同定領域の表面部に強制変位の拘束条件を割り当て、前記材料モデル同定領域の内部に強制変位の拘束条件を割り当てず、前記非同定領域の表面部及び内部に強制変位の拘束条件を割り当てる、請求項８記載の血管解析装置。
  10. 　前記力学モデルに基づく、血圧による内圧の時系列の空間分布、材料モデルパラメータ、及び血管形態指標を関連付けて記憶する関連付け記憶部をさらに備える、請求項１記載の血管解析装置。
  11. 　前記時系列の形態指標と前記時系列の形状変形指標とは確率分布で与えられる、請求項１記載の血管解析装置。
  12. 　前記構築部は、前記暫定的に構築された力学モデルを超音波診断装置からのプラーク指標に応じて複数の部分領域に区分し、前記複数の部分領域の各々に材料モデルのためのパラメータを割り当てる、請求項１記載の血管解析装置。
  13. 　前記血液流量指標の観測値は、前記時系列の医用画像に基づく血管内腔の容積の時間変化、あるいは、超音波ドプラ計測に基づく血流量の時間変化である、請求項１記載の血管解析装置。
  14. 　被検体の血管に関する時系列の医用画像のデータを発生する撮像機構と、
    　前記時系列の医用画像に含まれる血管領域に解析対象領域を設定し、前記解析対象領域に潜在変数同定領域を設定する設定部と、
    　前記時系列の医用画像を画像処理して前記解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出する算出部と、
    　前記時系列の形態指標と前記時系列の形状変形指標と前記時系列の医用画像とに基づいて前記解析対象領域の構造流体解析に関する力学モデルを暫定的に構築する構築部と、
    　前記暫定的に構築された力学モデルに基づく血管形態指標の予測値及び血液流量指標の予測値の少なくとも一方が予め計測された血管形態指標の観測値及び血液流量指標の観測値の少なくとも一方に整合するように前記潜在変数同定領域に関する潜在変数を同定する同定部と、
    　を具備する医用画像診断装置。
  15. 　時系列の医用画像に含まれる血管領域に解析対象領域を設定し、前記解析対象領域に潜在変数同定領域を設定し、
    　前記時系列の医用画像を画像処理して前記解析対象領域の時系列の形態指標と時系列の形状変形指標とを算出し、
    　前記時系列の形態指標と前記時系列の形状変形指標と前記時系列の医用画像とに基づいて前記解析対象領域の構造流体解析に関する力学モデルを暫定的に構築し、
    　前記暫定的に構築された力学モデルに基づく血管形態指標の予測値及び血液流量指標の予測値の少なくとも一方が予め計測された血管形態指標の観測値及び血液流量指標の観測値の少なくとも一方に整合するように前記潜在変数同定領域に関する潜在変数を同定する、
    　ことを具備する血管解析方法。
PCT/JP2013/080502 2012-11-19 2013-11-12 血管解析装置、医用画像診断装置、及び血管解析方法 WO2014077233A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
CN201380054742.6A CN104736061A (zh) 2012-11-19 2013-11-12 血管解析装置、医用图像诊断装置及血管解析方法
US14/716,433 US11357409B2 (en) 2012-11-19 2015-05-19 Blood vessel analysis apparatus, medical image diagnosis apparatus, and blood vessel analysis method

Applications Claiming Priority (2)

Application Number Priority Date Filing Date Title
JP2012253619A JP2014100249A (ja) 2012-11-19 2012-11-19 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム
JP2012-253619 2012-11-19

Related Child Applications (1)

Application Number Title Priority Date Filing Date
US14/716,433 Continuation US11357409B2 (en) 2012-11-19 2015-05-19 Blood vessel analysis apparatus, medical image diagnosis apparatus, and blood vessel analysis method

Publications (1)

Publication Number Publication Date
WO2014077233A1 true WO2014077233A1 (ja) 2014-05-22

Family

ID=50731144

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2013/080502 WO2014077233A1 (ja) 2012-11-19 2013-11-12 血管解析装置、医用画像診断装置、及び血管解析方法

Country Status (4)

Country Link
US (1) US11357409B2 (ja)
JP (1) JP2014100249A (ja)
CN (2) CN107582041B (ja)
WO (1) WO2014077233A1 (ja)

Cited By (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2017516535A (ja) * 2014-05-30 2017-06-22 ハートフロー， インコーポレイテッド 血流特性の報告のためのシステム及び方法
EP3220811A4 (en) * 2014-11-17 2018-07-25 David A. Borkholder Blood pressure and arterial compliance estimation from arterial segments
CN109478332A (zh) * 2016-05-16 2019-03-15 阿西斯特医疗系统有限公司 基于运动的图像分割系统和方法

Families Citing this family (45)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2013507227A (ja) 2009-10-12 2013-03-04 シリコンバレー メディカル インスツルメンツ インコーポレイテッド コレジスタ・イメージングのための血管内超音波システム
JP2014100249A (ja) 2012-11-19 2014-06-05 Toshiba Corp 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム
JP6091870B2 (ja) 2012-12-07 2017-03-08 東芝メディカルシステムズ株式会社 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム
EP2784748B1 (en) * 2013-03-28 2017-11-01 Expert Ymaging, SL A computer implemented method for assessing vascular networks from medical images and uses thereof
US9986938B2 (en) * 2014-02-25 2018-06-05 Medis Associated B.V. Method and device for determining a geometrical parameter of a blood vessel
NL2012459B1 (en) 2014-03-18 2016-01-08 Medis Ass B V Method and device for determining deviation in pressure in a blood vessel.
US9058692B1 (en) * 2014-04-16 2015-06-16 Heartflow, Inc. Systems and methods for image-based object modeling using multiple image acquisitions or reconstructions
US10636146B2 (en) * 2015-05-12 2020-04-28 Singapore Health Services Pte Ltd Medical image processing methods and systems
EP3310261A4 (en) 2015-06-19 2019-01-23 Neural Analytics, Inc. TRANSCRANIAL DOPPLER PROBE
US10909661B2 (en) 2015-10-08 2021-02-02 Acist Medical Systems, Inc. Systems and methods to reduce near-field artifacts
US10653393B2 (en) 2015-10-08 2020-05-19 Acist Medical Systems, Inc. Intravascular ultrasound imaging with frequency selective imaging methods and systems
US11369337B2 (en) * 2015-12-11 2022-06-28 Acist Medical Systems, Inc. Detection of disturbed blood flow
US11589836B2 (en) 2016-01-05 2023-02-28 Novasignal Corp. Systems and methods for detecting neurological conditions
JP2019500155A (ja) 2016-01-05 2019-01-10 ニューラル アナリティクス、インコーポレイテッド 一体型プローブ構造
WO2017120388A1 (en) 2016-01-05 2017-07-13 Neural Analytics, Inc. Systems and methods for determining clinical indications
US10278662B2 (en) 2016-02-05 2019-05-07 Toshiba Medical Systems Corporation Image processing apparatus and medical image diagnostic apparatus
JP2017176202A (ja) * 2016-03-28 2017-10-05 コニカミノルタ株式会社 動態解析システム
JP6378715B2 (ja) * 2016-04-21 2018-08-22 ゼネラル・エレクトリック・カンパニイ 血管検出装置、磁気共鳴イメージング装置、およびプログラム
JP6929689B2 (ja) * 2016-04-26 2021-09-01 キヤノンメディカルシステムズ株式会社 医用画像処理装置及び医用画像診断装置
US9947093B2 (en) 2016-05-03 2018-04-17 Konica Minolta, Inc. Dynamic analysis apparatus and dynamic analysis system
US10610184B2 (en) 2016-08-12 2020-04-07 Canon Medical Systems Corporation Medical-information processing apparatus and X-ray CT apparatus
DE202017106016U1 (de) 2016-10-04 2017-12-07 Toshiba Medical Systems Corporation Medizinische Informationsverarbeitungsvorrichtung, Röntgen-CT-Vorrichtung und computerlesbares Speichermedium mit einem Programm für ein medizinisches Informationsverarbeitungsverfahren
JP2018068631A (ja) * 2016-10-28 2018-05-10 キヤノン株式会社 放射線撮影システム、放射線表示方法
US10163209B2 (en) * 2016-11-23 2018-12-25 Toshiba Medical Systems Corporation Medical image processing apparatus, medical image processing method, and X-ray CT apparatus
US11202612B2 (en) * 2016-11-28 2021-12-21 Canon Medical Systems Corporation Medical image-processing apparatus, X-ray CT apparatus, and medical image-processing method performing fluid analysis to switch displayed color information
CN106780527B (zh) * 2016-11-29 2020-09-15 上海联影医疗科技有限公司 医学图像中血管进出口、边界条件获取方法及处理装置
US11517277B2 (en) * 2016-12-15 2022-12-06 Koninklijke Philips N.V. Visualizing collimation errors
WO2018133118A1 (zh) 2017-01-23 2018-07-26 上海联影医疗科技有限公司 血流状态分析系统及方法
EP3382583A1 (en) * 2017-03-29 2018-10-03 Koninklijke Philips N.V. Hemodynamic simulation of movement inducted vascular deformations
JP6847453B2 (ja) * 2017-05-30 2021-03-24 国立大学法人九州大学 地図情報生成方法、作動方法、及びプログラム
JP7073661B2 (ja) * 2017-09-27 2022-05-24 コニカミノルタ株式会社 動態解析装置及び動態解析システム
EP3691531A4 (en) * 2017-10-06 2021-05-26 Emory University METHODS AND SYSTEMS FOR DETERMINING HEMODYNAMIC INFORMATION FOR ONE OR MORE ARTERIAL SEGMENTS
JP7170459B2 (ja) * 2018-08-10 2022-11-14 株式会社東芝 血行検出装置、方法、及びプログラム
CN109620199B (zh) 2018-11-30 2021-03-16 博动医学影像科技（上海）有限公司 建立血管截面函数、血管压力差和血管应力的方法及装置
EP3685756A1 (en) * 2019-01-24 2020-07-29 Koninklijke Philips N.V. Methods and systems for investigating blood vessel characteristics
JP2019171102A (ja) * 2019-05-30 2019-10-10 コニカミノルタ株式会社 画像処理装置、画像処理方法、およびプログラム
US11024034B2 (en) 2019-07-02 2021-06-01 Acist Medical Systems, Inc. Image segmentation confidence determination
CN110458848B (zh) * 2019-07-05 2021-11-26 心医国际数字医疗系统(大连)有限公司 判断影像的冠脉分割中冠脉分支的方法
CN110559015B (zh) * 2019-08-26 2020-12-22 清华大学 血管生理参数测量方法、设备、计算机设备和存储介质
CN112151180B (zh) * 2019-12-05 2024-03-08 苏州润迈德医疗科技有限公司 具有狭窄病变的血管数学模型的合成方法和装置
US20220110530A2 (en) * 2019-12-09 2022-04-14 Nordsletten David Method and System for Estimating Pressure Difference in Turbulent Flow
CN113876304A (zh) * 2021-09-08 2022-01-04 深圳市中科微光医疗器械技术有限公司 一种基于oct图像和造影图像确定ffr的方法及装置
JP2023061477A (ja) * 2021-10-20 2023-05-02 富士通株式会社 プログラム、データ処理方法及びデータ処理装置
CN116705330B (zh) * 2023-07-31 2023-11-10 柏意慧心(杭州)网络科技有限公司 确定血管壁的弹性特征的方法、计算设备和介质
CN117316359B (zh) * 2023-09-22 2024-04-09 杭州威灿科技有限公司 血液检测过程跟踪方法、装置、设备及介质

Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2008514368A (ja) * 2004-09-30 2008-05-08 ザ リージェント オブ ザ ユニバーシティ オブ カリフォルニア 特性判定方法、記録装置、評価方法
WO2012028190A1 (en) * 2010-09-02 2012-03-08 Pie Medical Imaging Bv Method and apparatus for quantitative analysis of a tree of recursively splitting tubular organs

Family Cites Families (13)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
US20050159896A1 (en) * 2002-04-04 2005-07-21 Toshio Ishikawa Apparatus and method for analyzing data
US7117026B2 (en) * 2002-06-12 2006-10-03 Koninklijke Philips Electronics N.V. Physiological model based non-rigid image registration
WO2006039358A2 (en) 2004-09-30 2006-04-13 The Regents Of The University Of California Method for assessment of the structure-function characteristics of structures in a human or animal body
WO2007066249A2 (en) 2005-12-09 2007-06-14 Koninklijke Philips Electronics, N.V. Model-based flow analysis and visualization
RU2455689C2 (ru) * 2006-11-22 2012-07-10 Конинклейке Филипс Электроникс Н.В. Создание изображения на основе ограниченного набора данных
CN101680018B (zh) * 2007-01-10 2017-03-15 海莫希尔有限责任公司 体外血液动力学的内皮/平滑肌细胞共培养模型在鉴定血管疾病的新型治疗靶标中的应用
JP4439533B2 (ja) 2007-03-27 2010-03-24 株式会社東芝 負荷算定装置および負荷算定方法
EP2160135B1 (en) 2007-06-21 2012-08-15 Koninklijke Philips Electronics N.V. Adjusting acquisition protocols for dynamic medical imaging using dynamic models
US8675940B2 (en) 2009-10-27 2014-03-18 Siemens Aktiengesellschaft Generation of moving vascular models and blood flow analysis from moving vascular models and phase contrast MRI
US8157742B2 (en) 2010-08-12 2012-04-17 Heartflow, Inc. Method and system for patient-specific modeling of blood flow
US8315812B2 (en) 2010-08-12 2012-11-20 Heartflow, Inc. Method and system for patient-specific modeling of blood flow
JP2014100249A (ja) * 2012-11-19 2014-06-05 Toshiba Corp 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム
JP6091870B2 (ja) * 2012-12-07 2017-03-08 東芝メディカルシステムズ株式会社 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム

Patent Citations (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2008514368A (ja) * 2004-09-30 2008-05-08 ザ リージェント オブ ザ ユニバーシティ オブ カリフォルニア 特性判定方法、記録装置、評価方法
WO2012028190A1 (en) * 2010-09-02 2012-03-08 Pie Medical Imaging Bv Method and apparatus for quantitative analysis of a tree of recursively splitting tubular organs

Cited By (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2017516535A (ja) * 2014-05-30 2017-06-22 ハートフロー， インコーポレイテッド 血流特性の報告のためのシステム及び方法
JP2018140227A (ja) * 2014-05-30 2018-09-13 ハートフロー， インコーポレイテッド 血流特性の報告のためのシステム及び方法
EP3220811A4 (en) * 2014-11-17 2018-07-25 David A. Borkholder Blood pressure and arterial compliance estimation from arterial segments
US11622730B2 (en) 2014-11-17 2023-04-11 Rochester Institute Of Technology Pulse wave velocity, arterial compliance, and blood pressure
CN109478332A (zh) * 2016-05-16 2019-03-15 阿西斯特医疗系统有限公司 基于运动的图像分割系统和方法

Also Published As

Publication number Publication date
CN107582041B (zh) 2021-07-06
JP2014100249A (ja) 2014-06-05
US11357409B2 (en) 2022-06-14
CN104736061A (zh) 2015-06-24
CN107582041A (zh) 2018-01-16
US20150245776A1 (en) 2015-09-03

Similar Documents

Publication Publication Date Title
WO2014077233A1 (ja) 血管解析装置、医用画像診断装置、及び血管解析方法
JP6918912B2 (ja) 画像処理装置、画像処理方法、及びプログラム
JP6362853B2 (ja) 血管解析装置、および血管解析装置の作動方法
US11887305B1 (en) Systems and methods for identifying personalized vascular implants from patient-specific anatomic data
US10206587B2 (en) Image processing apparatus, image processing method, and storage medium
JP6091870B2 (ja) 血管解析装置、医用画像診断装置、血管解析方法、及び血管解析プログラム
JP6275797B2 (ja) 管状構造解析装置、管状構造解析方法及び管状構造解析プログラム
JP6362851B2 (ja) 血管解析装置、血管解析プログラム、及び血管解析装置の作動方法
KR102401542B1 (ko) 관동맥 병변들의 위치, 발병, 및/또는 변화를 예측하기 위한 시스템들 및 방법들
JP6818492B2 (ja) 画像処理装置、画像処理方法、及びプログラム
JP6425916B2 (ja) 処理装置、画像処理方法、プログラム、および、造形装置
WO2013180851A1 (en) A method for estimating pressure gradients and fractional flow reserve from computed tomography angiography: transluminal attenuation flow encoding
US20220151500A1 (en) Noninvasive quantitative flow mapping using a virtual catheter volume
Chen et al. A personalized pulmonary circulation model to non-invasively calculate fractional flow reserve for artery stenosis detection
JP6675458B2 (ja) 血管解析装置、血管解析方法及び血管解析プログラム
JP2015217113A (ja) 血管解析装置、医用画像診断装置、血管解析方法及び血管解析プログラム
Kolli et al. Improved Functional Assessment of Ischemic Severity Using 3D Printed Models
Yamaguchi Computational Mechanical Study of the Blood Flow for the Application to Clinical Cardiovascular Medicine
Graf et al. Clinical benefits of integrating cardiac and vascular models

Legal Events

Date Code Title Description
121 Ep: the epo has been informed by wipo that ep was designated in this application

Ref document number: 13855415

Country of ref document: EP

Kind code of ref document: A1

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 13855415

Country of ref document: EP

Kind code of ref document: A1