JP5643580B2 - Blood flow dynamic analysis device, blood flow dynamic analysis program, fluid analysis device, and fluid analysis program - Google Patents

Description

  Embodiments described herein relate generally to a blood flow dynamic analysis device, a blood flow dynamic analysis program, a fluid analysis device, and a fluid analysis program.

  Conventionally, X-ray CT (Computed Tomography) devices and MRI (Magnetic Resonance Imaging) devices administer contrast agents from veins, collect time-series image data, analyze them based on the images, and blood flow in tissues. Tests are being conducted to obtain information on dynamics. This is called a perfusion examination, and utilizes the fact that the degree of concentration of a contrast agent on an imaging section can be collected as a change in image density.

For example, in order to eliminate variations in pulmonary circulation and contrast agent administration, a time concentration curve (TCC) of an artery flowing into a tissue: _ai as an input function and a decontour with a measured tissue time concentration curve C _i Deconvolution (Deconvolution) is performed to obtain a tissue-specific impulse response function (transfer function) R _i (residue function). From the obtained impulse response function, blood, which is an index that quantitatively represents hemodynamics Parameters such as flow rate (CBF: Cerebral Blood Flow in the case of brain), average transit time: MTT (Mean Transit Time), and blood volume (CBV: Cerebral Blood Volume in the case of brain) are calculated. This analysis method is called a deconvolution method.

  The time density curve is a curve representing a time change of an amount proportional to the contrast agent concentration. A time change curve of an amount proportional to the contrast agent concentration obtained from the contrast image data collected by the X-ray CT apparatus is also referred to as TDC (Time Density Curve), but is expressed as TCC here.

  For example, in Non-Patent Document 1, a standard SVD (Singular Value Decomposition) method or a block circulant SVD method is applied to deconvolution.

Here, the standard SVD method will be described. First, in the collected image data, a change in the pixel value of the pixel j (position r _j ) is measured to calculate a tissue temporal density curve C _i (r _j ). Also, an arterial time concentration curve A _i is calculated from the calculated tissue time concentration curve C _i as shown in the following equation (1). In Expression (1), i represents a phase number, and Ω _A represents a set of pixels corresponding to an artery.

Assuming that the tissue time concentration curve C _i is a deconvolution (convolution integration) of the arterial time concentration curve A _i and the impulse response function X, it can be expressed by the following system equation. In the formula (2), T _R is the time interval of image capturing (repetition time TR: repetition time) represents, R _j represents the impulse response function.

For example, the arterial time concentration curve a _i and the tissue time concentration curve C _i are discretized in the time axis direction, that is, the phase range of i _bs ≦ i ≦ L−1 is divided into L equal intervals, Assuming that data of L phases are used for analysis, the above equation (2) becomes C = AX in matrix form and can be expressed as shown in the following equation (3). In Equation (3), C is an L × 1 column vector, A is an L × L coefficient matrix, and X is an L × 1 column vector.

そして、上記の式（３）の近似解を、特異値分解により求める。すなわち、行列Ａを次式（５）に示すように３つの行列に分解することができる。
［数５］
Ａ＝ＵＤＶ^Ｔ ・・・（５） In the above equation (3), i _bs represents the first phase of the phase range to be analyzed, and the phase in which the contrast density of the artery and tissue starts to rise or a phase before that is used. Each element value a _{i, j} of the coefficient matrix is expressed by the following equation (4). In equation (4), i = 0, 1,... L-1, and j = 0, 1,.

And the approximate solution of said Formula (3) is calculated | required by singular value decomposition. That is, the matrix A can be decomposed into three matrices as shown in the following equation (5).
[Equation 5]
A = UDV ^T (5)

  In the above equation (5), U is an m-by-k orthonormal matrix, D is a k-by-k diagonal matrix, and V is an n-by-k orthonormal matrix. When used in the above equation (3), m = n = k = L. The diagonal matrix is a matrix in which all components other than the diagonal component are 0, and the orthonormal matrix is a diagonal matrix in which the diagonal component is 1.

_{Here, U = (u 1, u} 2, ... u k), V = (v 1, v 2, ... v k), D = diag (w 0, w 1, ... w k ), W _k is rearranged in descending order, and u _k and v _k are also replaced so as to correspond to the order. Using the result of the singular value decomposition of the matrix A, the minimum norm estimated value that is the least square solution of X can be calculated by the following equation (6).

上記の式（７）で算出されたインパルス応答関数Ｒ_ｉから、組織血流量、組織血液量、および平均通過時間を算出することができる。 In the above equation (6), W is a diagonal matrix composed of reciprocals of diagonal elements (singular values) of D, and W = diag (1 / w _o , 1 / w ₁ ,... 1 w _G−1 , 0, 0,. Here, for a singular value that is less than P _svd times the maximum singular value, 0 is used as a diagonal element instead of the reciprocal number. Accordingly, only G diagonal elements have non-zero values. P _svd is a parameter that influences the strength of regularization. The larger the value, the stronger regularization is performed.
When the norm estimated value is obtained, the impulse response function R _i can be calculated as shown in the following equation (7).

From the impulse response function R _i calculated by the above equation (7), the tissue blood volume, the tissue blood volume, and the average transit time can be calculated.

  By the way, the blood flow flows into the capillaries of the tissue, but the time for the contrast medium to reach the tissue differs depending on the location of the tissue, and the time density curve of the artery obtained from the set arterial region and the time density of the tissue of each pixel Strictly speaking, the time axis does not coincide. Therefore, in extreme cases, before the contrast density of the arteries rises (before the arterial time density curve rises), the tissue contrast density rises (the tissue time density curve rises first). It becomes.

  In the standard SVD method described above, it is assumed that there is no such time-direction shift (delay). However, there is always a time lag between the arterial time concentration curve and the tissue time concentration curve, so the deconvolution result varies depending on the delay time value, and the calculated tissue blood flow value has an effect. You will receive it.

Therefore, in the block circulant SVD method, a system matrix that is a block circulant matrix as shown in the following equation (8) is used in order to reduce the dependency on the deviation in the time direction.

  In this block circulant SVD method, the basic processing is the same as in the standard SVD method, and the matrix format C = AX, but A is a 2L × 2L coefficient matrix. Therefore, 0 is added to the back of each element value of the coefficient matrix A constituting the arterial time density curve, and 0 is added to the back of each element value of the column vector X constituting the impulse response function. 0 is added to the back of each element value of the column vector C constituting the curve, and the number of elements is doubled to the measurement data to be used. In particular, the coefficient matrix A is rearranged so that the rows and columns are both doubled or expanded, and the elements of the column vectors are periodic.

  In this way, by configuring the system matrix as shown in the above equation (8), when there is a time direction deviation, the impulse response function is shifted in the time axis direction so that the shape thereof maintains a similar shape. become. Therefore, by shifting the impulse response function in the time axis direction, it is possible to calculate the tissue blood flow volume, the tissue blood volume, and the average transit time with less dependency on the shift in the time direction.

  However, in the system matrix, the coefficient matrix A constituting the arterial time density curve is 2 × 2 times, and the column vector constituting the tissue time density curve is twice, so that the analysis processing time is 8 times or more. It increases and is not practical.

Therefore, in Non-Patent Document 2, the system matrix of the standard SVD method is applied without increasing the analysis processing time. That is, in the method proposed in Non-Patent Document 2, the above equation (4) is used as the system matrix, and the elements of the time density curve of the tissue from the first phase to the offset position are included in the known vector as noise (signal (Strain) is excluded from the calculation. Here, in the above equation (3), since the leading number phase from the vector C to the offset position is removed, the same number of zero elements is added behind.
As described above, in the technique of Non-Patent Document 2, the dependency on the shift in the time direction is reduced without expanding the system matrix.

Ona Wu et al, "Tracer Arrival Timing-Insensitive Technique for Estimating Flow in MR Perfusion-Weighted Imaging Using Singular Value Decomposition With a Block-Circulant Deconvolution Matrix", Magnetic Resonance in Medicine 50: 164-174 (2003) M.R.Smith, H.Lu, S.Trochet, and R.Frayne "Removing the Effect of SVD Algorithmic Artifacts Present in Quantitative MR Perfusion Studies", Magnetic Resonance in Medicine 51: 631-634 (2004)

  In the technique of Non-Patent Document 1, the standard SVD method does not consider the deviation in the time direction, so the deconvolution result is not accurate. On the other hand, in the block circulant SVD method, although the shift in the time direction is taken into account, the analysis processing time is enormous, and there are problems that are not practical.

  In particular, in the operation of an X-ray CT apparatus and an MRI apparatus, it is usually necessary to examine a large number of patients. In addition, an analysis process must be performed immediately after the image is taken, and the engineer must check the result of the process to check whether or not the shooting is incomplete. Therefore, if the analysis processing time becomes enormous, the number of tests that can be performed decreases, the waiting time increases for the patient, and the acquisition of the patient decreases for the hospital, resulting in both disadvantages.

  Therefore, in the technique of Non-Patent Document 2, processing is performed without increasing the analysis processing time. From the calculation, the element of the time concentration curve of the tissue from the first phase to the offset position is regarded as noise (signal distortion). Since it excludes, the new problem that a result will depend on an offset position will arise. However, it is difficult to correctly identify the offset position, and as a result, there is a problem that information relating to blood flow dynamics cannot be obtained stably.

  The present invention has been made in view of such a situation, and an object of the present invention is to provide a blood flow dynamic analysis device and a blood flow dynamic that can stably calculate information related to fluid dynamics without increasing processing time. An analysis program, a fluid analysis apparatus, and a fluid analysis program are provided.

  A blood flow dynamic analysis apparatus according to an embodiment of the present invention includes a tissue time concentration curve calculation unit, an arterial time concentration curve calculation unit, a first configuration unit, a second configuration unit, a deconvolution integration unit, and blood flow information. Calculation means is provided. The tissue time density curve calculating means measures a change in the pixel value of each pixel between the plurality of time phases based on the contrast image data of the subject in the plurality of time phases obtained by the medical imaging apparatus. A time concentration curve representing a time change of the amount corresponding to the concentration of the contrast agent in is calculated. The arterial time concentration curve calculating means sets an arterial region from the time concentration curve corresponding to the tissue, and calculates a time concentration curve corresponding to the artery. The first component means discretizes the time concentration curve corresponding to the tissue with reference to the first time, and includes a plurality of elements including values corresponding to the concentration of the contrast agent in the tissue at different times Are arranged in a one-dimensional manner. The second component means discretizes the time concentration curve corresponding to the artery with reference to a second time behind the first time, and converts the time concentration curve to the concentration of the contrast agent in the artery at a different time. A second set in which a plurality of elements including corresponding values are two-dimensionally arranged is configured. The deconvolution integration means performs deconvolution integration of the first set and the second set to calculate a transfer function of the tissue. The blood flow information calculation means calculates information related to blood flow dynamics based on the transfer function.

  In addition, a blood flow dynamic analysis program according to an embodiment of the present invention includes a computer, a tissue time concentration curve calculation unit, an arterial time concentration curve calculation unit, a first configuration unit, a second configuration unit, It functions as a convolution integrator and a blood flow information calculator.

In addition, the fluid analysis apparatus according to the embodiment of the present invention includes a concentration change calculation unit and a fluid information calculation unit. The density change calculation means calculates the first change between the plurality of frames in an amount corresponding to the density of the pixel value in the pixel corresponding to the tissue based on the image data of the plurality of frames of the subject, while the fluid The second change between the plurality of frames is calculated in an amount corresponding to the density of the pixel value corresponding to. The fluid information calculation means creates the first plurality of elements to be analyzed based on the element corresponding to the first point of the first change based on the first change. , Based on the second change, the second plurality of the number to be analyzed with the element corresponding to the second point after the first point of the second change as the first element An element is created, and information on fluid dynamics is calculated with a deconvolution integration process based on the first plurality of elements and the second plurality of elements.
The fluid analysis program according to the embodiment of the present invention causes a computer to function as the above-described concentration change calculation unit and fluid information calculation unit.

The figure which shows the structural example of the blood-flow dynamics analysis apparatus which concerns on the 1st Embodiment of this invention. The block diagram which shows an example of the function of CPU which performs a local blood flow analysis application. The figure explaining a peak phase. The flowchart explaining the calculation process of the information regarding the blood flow dynamics of a structure | tissue. The figure which shows typically the system matrix which concerns on 1st Embodiment. The functional block diagram which shows the structural example of the fluid analysis apparatus which concerns on the 2nd Embodiment of this invention. The figure which shows the 1st example of the imaging conditions for non-contrast fluid imaging performed in the imaging part of the MRI apparatus shown in FIG. The figure which shows the 2nd example of the imaging conditions for non-contrast fluid imaging performed in the imaging part of the MRI apparatus shown in FIG. The figure which shows the 3rd example of the imaging conditions for non-contrast fluid imaging performed in the imaging part of the MRI apparatus shown in FIG. The figure explaining the creation method of the element used as the analysis object in the system matrix preparation part 3d shown in FIG. The flowchart which shows the flow of the data processing in the fluid analysis apparatus shown in FIG.

  Hereinafter, embodiments of the present invention will be described in detail with reference to the drawings.

(First embodiment)
FIG. 1 is a diagram illustrating a configuration example of a blood flow dynamic analysis apparatus 1 according to the first embodiment of the present invention.

  In the blood flow dynamic analysis device 1, a CPU (Central Processing Unit) 1a, a ROM (Read Only Memory) 1b, a RAM (Random Access Memory) 1c, and an input / output interface 1d are connected via a bus 1e. A storage unit 1f, a communication unit 1g, a display unit 1h, an input unit 1i, and a removable medium 1j are connected to the input / output interface 1d.

  The CPU 1a reads out and executes a boot program for starting the blood flow dynamic analysis device 1 from the ROM 1b based on an input signal from the input unit 1i, and reads out various operating systems stored in the storage unit 1f. The CPU 1a performs various controls based on an input signal from the input unit 1i, reads a program and data stored in the ROM 1b and the storage unit 1f and loads them into the RAM 1c, or a program read from the RAM 1c. Based on the command, a series of processes such as data calculation or processing is executed.

  The storage unit 1f includes a semiconductor memory, a magnetic disk, and the like, and stores programs and data executed by the CPU 1a. For example, an application for analyzing information related to local blood flow dynamics is prepared in the storage unit 1f as a program executed by the CPU 1a. Hereinafter, this application for analyzing information related to local blood flow dynamics is referred to as a local blood flow analysis application.

  The communication unit 1g is configured by a LAN (Local Area Network) card, a modem, or the like, and enables the blood flow dynamic analysis device 1 to be connected to a communication medium such as a LAN or the Internet. That is, the communication unit 1g transmits data received from the communication medium to the CPU 1a via the input / output interface 1d and the bus 1e, and transmits data received from the CPU 1a via the bus 1e and the input / output interface 1d to the communication medium. To do.

  The display unit 1h is, for example, a liquid crystal display. Based on signals received from the CPU 1a via the bus 1e and the input / output interface 1d, the processing result of the CPU 1a is displayed.

  The input unit 1i is configured by an input device such as a keyboard and a mouse through which an operator of the blood flow dynamic analysis apparatus 1 inputs various operations. An input signal is generated based on the operation of the operator and transmitted to the CPU 1a via the input / output interface 1d and the bus 1e.

  The removable medium 1j is, for example, an optical disk or a flexible disk. Data read by a disk drive (not shown) is transmitted to the CPU 1a via the input / output interface 1d and the bus 1e, and data received from the CPU 1a via the bus 1e and the input / output interface 1d is written by the disk drive. .

  The blood flow dynamics analysis apparatus 1 in the first embodiment analyzes blood flow dynamics from image data collected by dynamics imaging using a medical image imaging apparatus (so-called medical modality) such as an X-ray CT apparatus or an MRI apparatus. Device. Therefore, the blood flow dynamic analysis device 1 only needs to be in an environment where such image data can be obtained, and may be configured integrally with the medical image capturing device, or may be configured separately from the medical image capturing device. May be. When configured separately, the collected image data is sent to the blood flow dynamics analysis device 1 from a medical image photographing device (not shown) via the removable medium 1j or via the communication unit 1g.

  In addition, the blood flow dynamic analysis device 1 can transfer the blood flow dynamic analysis result to the PACS (Picture Archiving and Communication System) via the communication unit 1g and store it there. As a result, it is possible to easily search or browse the blood flow dynamic analysis result of specific image data.

Medical imaging apparatuses include, in addition to X-ray CT apparatuses and MRI apparatuses, ultrasonic diagnostic apparatuses and PET (Positron Emission Tomography) apparatuses.
FIG. 2 is a block diagram illustrating an example of the function of the CPU 1a that executes the local blood flow analysis application.

  As shown in FIG. 2, the CPU 1a that executes the local blood flow analysis application includes a time-series data acquisition unit 2a, a contrast density curve calculation unit 2b, an arterial time density curve calculation unit 2c, a system matrix creation unit 2d, and a deconvolution. It functions as part 2e.

The time-series data acquisition unit 2a acquires time-series imaging data (multi-time phase image data) generated by a medical imaging apparatus for a predetermined period of time for a subject into whom a contrast medium is injected intravenously (bolus injection). Then, the acquired time-series imaging data is output to the contrast density curve calculation unit 2b.
The contrast density curve calculation unit 2b receives time-series imaging data output from the time-series data acquisition unit 2a, measures changes in pixel values of each pixel, and calculates a contrast density curve.

  In other words, when a contrast medium is injected intravenously into a subject in a relatively short time, the contrast medium passes through the veins and returns to the heart, circulates through the lungs, returns to the heart, flows into the arteries, and reaches the tissues of each organ. Go back through the vein and back to the heart. Here, for example, when the medical image capturing apparatus captures a CT image or MRI image of the brain region in time series, it is possible to measure the change in the pixel value of each pixel according to the change in the contrast agent concentration. . The contrast density curve calculation unit 2b can calculate a contrast density curve based on the measured change in the pixel value of each pixel.

For example, in the case of a CT image, an amount proportional to the contrast agent concentration (TDC) is obtained by subtracting the pixel value S ₀ before arrival of the contrast agent from the pixel value S _i of each time phase as shown in the following equation (9). : Time Density Curve). In Expression (9), i represents a shooting time phase.
[Equation 9]
C _i = S _i −S ₀ (9)

Further, for example, in the case of an MRI image, as shown in the following equation (10), the contrast agent concentration is calculated from the logarithm (ΔR2 ^* ) of the ratio between the pixel value S _i of each time phase and the pixel value S ₀ before arrival of the contrast agent. A proportional amount (TCC: Time Concentration Curve) is calculated. In Expression (10), i represents the photographing time phase, and TE represents the echo time.
[Equation 10]
ΔR2 ^* = − log (S _i / S ₀ ) / TE (10)

Contrast concentration curve calculation unit 2b, the contrast concentration curve calculated by the formula (9) or formula (10), the time-density curve C _i of the tissue, the calculated value of the arterial time-density curve calculation unit 2c, and, The data is output to the system matrix creation unit 2d.

  The contrast density curve calculation unit 2b measures the change in the pixel value of each pixel between a plurality of time phases based on the contrast image data of the subject in a plurality of time phases, and an amount corresponding to the contrast agent concentration in the tissue It functions as a tissue time concentration curve calculating means for calculating a time concentration curve representing a time change of the tissue.

The arterial time density curve calculation unit 2c receives the tissue time density curve C _i output from the contrast density curve calculation unit 2b, and from the curve, the arterial of the voxel serving as an analysis unit of the CT image or the MRI image is calculated. A time concentration curve is calculated. In other words, the size of the voxel that is a unit of analysis for CT images and MRI images is about 0.5 mm to 3 mm, and if it is the main artery in the tissue, there are voxels that are almost included in the blood vessels. To do.

The arterial time concentration curve calculation unit 2c sets such an arterial region by specifying such voxels automatically or by designation by the user, and sets a time concentration curve (here, the i-th imaging) in the arterial region. It functions as an arterial time concentration curve calculating means for calculating the concentration value a _i ) at the time. When a plurality of pixels are specified as an arterial region, for example, the average value is set to a _i .

Alternatively, the arterial time concentration curve calculation unit 2c may convert a function of a constant multiple of the gamma probability density function according to the following equation (11) (a function proportional to the gamma probability density function) to the curve of the arterial time concentration curve as necessary. The model is applied as a model (hereinafter referred to as “arterial curve model”), and curve fitting (fitting to the arterial curve model) of the calculated arterial time concentration curve a _i is performed. Normally, curve fitting is often performed for MRI images, but curve fitting is often not performed for CT images.
[Equation 11]
φ (t; t ₀ , α, β, K) = K (tt ₀ ) ^α exp {-β ^-1 (tt ₀ )} (11)

  The gamma probability density function mentioned here represents the probability density of the time required to pass through all the multiple barriers α that can pass at a rate of once per time β. Since the arterial function can be thought of as a function distributed through the lung barrier, this gamma probability density function is applied as an arterial curve model to fit (approximate) the arterial curve model. .

The arterial time concentration curve a (t) is expressed by the following equation (12) by curve fitting using a function of a constant multiple of the gamma probability density function expressed by the above equation (11). In Expression (12), t represents time, and t ₀ , α, β, and K all represent parameters. α is the number of barriers, β is the time required to pass through all of the barriers α, and K is a linear parameter.
[Equation 12]
a (t) = φ (t; t ₀ , α, β, K) (12)

動脈の時間濃度曲線算出部２ｃは、カーブフィッティングした動脈の時間濃度曲線ａ_ｉの算出値を、システム行列作成部２ｄに出力する。 The arterial time concentration curve calculation unit 2c uses all the four parameters t ₀ , α, β, and K in the above equation (12) as variables, and the arterial curve model and the artery concentration value shown in the following equation (13): Calculate t ₀ , α, β, and K that minimize the square error with a _i . Actually, in many cases, a method of calculating α, β, and K using a specific time before the contrast medium reaches the target organ is used as t ₀ .

The arterial time concentration curve calculation unit 2c outputs the calculated value of the curve-fitted arterial time concentration curve a _i to the system matrix creation unit 2d.

The system matrix creation unit 2d receives the value of the tissue time density curve C _i output from the contrast density curve calculation unit 2b and the arterial time density curve a _i output from the arterial time density curve calculation unit 2c. Enter a value for.

The system matrix creation unit 2d determines the first phase i _bs and phase range i _bs ≦ i ≦ L−1 to be analyzed from the input time concentration curve C _i of the tissue and the time concentration curve a _{i of the} artery. As the first phase i _bs to be analyzed, a phase in which the contrast density of the artery and tissue starts to rise or a phase before that is used. A phase is a time determined by a sampling interval in analysis, and a phase range is a so-called analysis range.

The system matrix creation unit 2d discretizes the arterial time concentration curve a _i in the determined phase range in the time axis direction, that is, divides the phase range of i _bs ≦ i ≦ L−1 into L equal intervals. The L pieces of data thus calculated are calculated as shown in the following equation (14).
[Formula 14]
A _i = a (iT _R ) (14)

In the above formula (14), T _R is represents the time interval of image capturing (repetition time) is not limited to this, a dynamic time phases, or the actual analysis time phase obtained from the dynamic time phase But you can. Further, the measurement value a _i may be used as A _i as it is.

System matrix creating unit 2d was obtained by discretizing the time-density curve C _i of the tissue in the time axis direction, the set consisting of density values of the tissue at different times, in front of the time axis direction, the H-number Add an element.

As a result, H elements are added in front of the leading phase i _bs, and a plurality of tissue concentration values at different times are obtained with reference to the collection start time or a predetermined time (the Nth phase from the collection start time). A first set of elements is arranged in one dimension (time axis direction), and a (H + L) × 1 column vector is configured.

  The system matrix creation unit 2d discretizes the time density curve corresponding to the tissue with reference to the first time (the collection start time or the Nth phase from the collection start time), and the contrast medium in the tissue at different times. It functions as a first component that constitutes a first set in which a plurality of elements including values corresponding to the density are arranged one-dimensionally.

In addition, the system matrix creation unit 2d generates a set of arterial concentration values at different times, obtained by discretizing the arterial time concentration curve A _i in the time axis direction, at the rear of the time axis direction. Add elements. As a result, the second set of a plurality of elements including the arterial concentration values at different times is arranged in one dimension (in the direction of the time axis) to form a (L + H) × 1 column vector. Then, the system matrix creation unit 2d constructs a (L + H) × (L + H) coefficient matrix in which the column vectors of the arterial density values arranged one-dimensionally are further arranged two-dimensionally (in the spatial direction).

As a result, the second set of a plurality of elements including the arterial density values at different times is arranged in two dimensions (time axis direction and spatial direction) with reference to the first phase i _bs to be analyzed. A coefficient matrix of (L + H) × (L + H) is constructed.

  The system matrix creation unit 2d generates a time concentration curve corresponding to the artery at a second time (the time concentration curve of the artery increases after the first time (collection start time or N phase from the collection start time)). Functioning as a second configuration means that constitutes a second set in which a plurality of elements including a value corresponding to the contrast agent concentration in the artery at different times are arranged two-dimensionally. To do.

Since the deconvolution of the arterial time concentration curve and the impulse response function is a tissue time concentration curve, the system matrix creation unit 2d has a plurality of elements including the tissue concentration values at different times configured as described above. C = AX as shown in the following equation (15), using a column vector C consisting of: a coefficient matrix A consisting of a plurality of elements including arterial density values at different times, and an unknown number X (impulse response function). Create a system matrix. In Expression (15), C is an (H + L) × 1 vector, A is an (L + H) × (L + H) matrix, and X is an (L + H) × 1 vector. i _bs represents the first phase to be analyzed.

  In the above equation (15), H represents the number of elements added to the L × L coefficient matrix. For example, the number obtained by dividing the range from the imaging (collection) start time to the first phase by the phase interval, or more or less may be used. If the number H of additional elements is reduced, the size of the matrix to be analyzed is reduced and the processing time is shortened. However, if the number is too small, the dependency on the deviation in the time direction increases. Therefore, if the number of values obtained by dividing the range from the photographing start time to the first phase by the phase interval is used, the value is sufficient to reduce the dependency on the deviation in the time direction.

Among the plurality of elements obtained by discretizing the time density curve of the tissue in the time axis direction, an element of 0 is input (addition) to the additional element value P _i (i = 0, 1,... H−1). ), but, partially or entirely as C _i, the measured tissue may be used the values of the time-density curve. In this case, actually, it is necessary to capture image data of L or more phases.

An element value a _{i, j obtained} by discretizing the arterial time density curve in the time axis direction is expressed by the following equation (16). In equation (16), i = 0, 1,... L-1, and j = 0, 1,.

システム行列作成部２ｄは、作成した上記の式（１５）をデコンボリューション部２ｅに出力する。 In the above equation (16), an element of 0 is input (added) to the additional element value q _{i, j} , but a part or the whole is measured as shown in the following equation (17). The value of the arterial time density curve may be used. In this case, actually, it is necessary to capture image data of L or more phases.

The system matrix creation unit 2d outputs the created formula (15) to the deconvolution unit 2e.

The deconvolution unit 2e functions as a deconvolution integration unit that receives the coefficient matrix output from the system matrix creation unit 2d and performs deconvolution using a method such as the SVD method or the Fourier transform method. Deconvolution unit 2e, from the unknowns X obtained by deconvolution, calculates the impulse response function R _i as shown in equation (7).

The deconvolution unit 2e searches the existence range of the peak phase (maximum value of the impulse response function) from the impulse response function R _i calculated by the above equation (7), and the peak range (width) from the searched peak phase. Is calculated.
FIG. 3 is a diagram for explaining the peak phase.

As shown in FIG. 3, the peak of the impulse response function R _i calculated in the first embodiment is a phase that is approximately H backwards compared to the conventional standard SVD method. Therefore, the deconvolution unit 2e must take into account that the peak phase is shifted backward, and therefore it is necessary to determine the peak range in advance before calculating information on local blood flow dynamics.

  As a calculation method of the peak range, for example, the peak phase is searched in the time axis direction from the first phase, and the peak range is set up to the phase having a minimum negative value before and after the searched peak phase. In addition, you may make it make it a peak range not to the negative minimum value but to the phase which takes the value of 0.

The deconvolution unit 2e, based on the searched peak phase and peak range, is information on local blood flow dynamics, such as mean transit time MTT (Mean Transit Time), local cerebral blood flow CBF (Cerebral Blood Flow), and local brain A blood volume CBV (Cerebral Blood Volume) is calculated as shown in the following equation (18). That is, the local cerebral blood flow rate CBF is obtained from the peak height of the impulse response curve, the local cerebral blood volume CBV is obtained from the peak range area of the impulse response curve, and the average transit time MTT is obtained from the peak range of the impulse response curve.
[Equation 18]
CBF = MAX (R _i )
CBV = SUM (R _i )
MTT = SUM (R _i ) / MAX (R _i ) (18)

  The deconvolution unit 2e also functions as blood flow information calculation means for calculating information related to blood flow dynamics based on a tissue transfer function (impulse response function).

The deconvolution unit 2e corrects the calculated average transit time MTT, local cerebral blood flow CBF, and local cerebral blood volume CBV as shown in the following equation (19).
[Equation 19]
CBV = (100h / ρ) T _R SUM (R _i )
CBF (r) = 60 (100h / ρ) MAX (R _i )
MTT2 (r) = T _R {SUM (R _i ) / MAX (R _i )} (19)

In the above equation (19), ρ represents the specific gravity of the brain, and h represents the hematocrit value (ratio of red blood cells contained in a certain amount of blood). The local cerebral blood flow CBF may be calculated (corrected) using the relationship of the following equation (20).
[Equation 20]
CBF = CBV / MTT2 (20)

Next, with reference to the flowchart of FIG. 4, the calculation process of the information regarding the blood flow dynamics of the tissue in the first embodiment will be described. In this process, in particular, a brain tissue is taken as an example, and a process for calculating information related to blood flow dynamics will be described.
In step S1, the time series data acquisition unit 2a acquires time series imaging data generated by imaging for a predetermined time by the medical imaging apparatus.

In step S2, the contrast density curve calculation unit 2b calculates a contrast density curve by measuring a change in the pixel value of each pixel from the time-series imaging data acquired in step S1. The contrast concentration curve, and the time-density curve C _i of the organization.

In step S3, the time-density curve calculation unit 2c of the artery, the time-density curve C _i of the tissue that has been calculated in the processing in step S2, sets the arterial region to be analyzed unit of CT image and MRI image, the artery region The time density curve a _{i at} is calculated. At this time, when the acquired imaging data is an MRI image, the arterial time concentration curve calculation unit 2c applies a function of a constant multiple of the gamma probability density function as an arterial curve model, and calculates the calculated arterial time concentration curve a _i . Perform curve fitting.
In step S4, the system matrix creation unit 2d determines the first phase i _bs and the phase range i _bs ≦ i ≦ L−1 to be analyzed.

In step S5, the system matrix creation unit 2d determines the tissue time in the phase range i _bs ≦ i ≦ L−1 determined in the process of step S4 out of the tissue time concentration curve C _i calculated in the process of step S2. A density vector is discretized in the time axis direction, and a set of tissue density values at different times obtained by discretization is arranged in one dimension to form a column vector.

In step S5, the system matrix creation unit 2d determines the artery in the phase range i _bs ≦ i ≦ L−1 determined in the process of step S4 from the time density curve a _i of the artery calculated in the process of step S3. The time density curve is discretized in the time axis direction, and a set of arterial density values at different times obtained by discretization is two-dimensionally arranged to form a coefficient matrix.

In step S6, the system matrix creation unit 2d adds H elements ahead in the time axis direction of a column vector (set) composed of tissue concentration values at different times configured in the process of step S5. As a result, the first phase i _{bs to} be analyzed is shifted by H in the time axis direction, and the first phase consisting of a plurality of elements including tissue concentration values at different times with reference to the collection start time or the predetermined time. The set is arranged one-dimensionally, and a column vector C of (H + L) × 1 is configured.

  In step S6, the system matrix creation unit 2d adds H elements to the rear of the coefficient matrix (set) composed of arterial density values at different times configured in step S5 in the time axis direction. To do. As a result, the second set of a plurality of elements including the concentration values of the artery at different times is arranged two-dimensionally with reference to the phase in which the contrast density of the artery starts to increase or the phase before that. , (L + H) × (L + H) coefficient matrix A is constructed.

In step S7, the system matrix creation unit 2d obtains a column vector C composed of a plurality of elements including tissue concentration values at different times and a plurality of arterial concentration values at different times obtained by the process of step S6. The system matrix shown in the above equation (15) is created using the coefficient matrix A composed of the following elements and the unknown number X (impulse response function).
FIG. 5 is a diagram schematically illustrating the system matrix of Expression (15) according to the first embodiment.

In FIG. 5, C schematically represents a column vector of (H + L) × 1 obtained by discretizing the tissue time concentration curve C _i in the time axis direction, and A is the arterial time concentration curve a _i. Is a schematic representation of a coefficient matrix of (L + H) × (L + H) discretized in the time axis direction, and X schematically represents an unknown (impulse response function).

As shown in the figure, in the column vector C constituting the time density curve of the tissue, H elements p _i are added in front, and in the coefficient matrix A constituting the time density curve of the artery, H elements are added behind. Elements q _{i, j} are added. As the number H of additional elements, it is preferable to use the value of the number obtained by dividing the range from the imaging start time to the first phase i _bs by the phase interval. The added element value P _i (i = 0, 1,... H−1) and the element value q _{i, j} are input (added) with an element of 0. The measured value may be used.

  By constructing such a system matrix, it is not necessary to double the coefficient matrix as in the conventional block circulant SVD method, and high-speed deconvolution processing is possible.

Returning to the description of FIG. In step S8, the deconvolution unit 2e performs deconvolution of the system matrix created in the process of step S7 using a method such as the SVD method or the Fourier transform method to obtain the vector X, and the above equation (7) The impulse response function R _i is calculated as shown in FIG.

In step S9, the deconvolution unit 2e searches for the peak phase from the impulse response function R _i calculated in the process of step S8, and calculates the peak range from the searched peak phase.

  In step S10, the deconvolution unit 2e, based on the peak phase and peak range calculated in step S9, information on local blood flow dynamics, that is, the average transit time MTT, the local cerebral blood flow CBF, and the local brain The blood volume CBV is calculated according to the above equation (18).

  In step S11, the deconvolution unit 2e corrects the average transit time MTT, the local cerebral blood flow CBF, and the local cerebral blood volume CBV calculated in the process of step S10 as shown in the above equation (19). .

  As described above, the blood flow dynamic analysis device 1 adds a predetermined number of elements in front of the column vector constituting the time concentration curve of the tissue to be analyzed, and behind the coefficient matrix constituting the arterial time concentration curve. By adding a predetermined number of elements to the first phase of the tissue time concentration curve and the reference position of the first phase of the arterial time concentration curve, the processing time is not increased and the shift in the time direction is prevented. It is possible to stably calculate the tissue blood flow volume, the tissue blood volume, and the average passage time, which are information related to the blood flow dynamics of the blood vessels, by reducing the dependency.

  As a result, the hospital side can carry out examinations of a large number of patients, and the examination waiting time of patients can be shortened. In addition, since the analysis result can be obtained stably, the situation such as re-examination (re-analysis) can be reduced, and the burden on the patient can be reduced.

  In addition, the blood flow dynamic analysis device 1 can display a corresponding image on the display unit 1h based on the calculated average transit time MTT, local cerebral blood volume CBV, and local cerebral blood flow CBF. It becomes possible to easily grasp the state of blood flow dynamics.

  In the above description, the brain tissue is described as an example of the target tissue. However, in the first embodiment, the target tissue is not limited to the brain tissue, and other tissues such as the liver, heart, or lung are targeted. Of course it is possible.

  In the above, correction of the hematocrit value and correction by the specific gravity of the brain are performed for the brain tissue. Therefore, for example, in the case where the tissue such as the liver or lung moves due to respiration or the heart is the target tissue, it is necessary to correct the motion of the tissue. In addition to the motion correction, it is necessary to appropriately correct the hematocrit value and other corrections.

(Second Embodiment)
FIG. 6 is a functional block diagram showing a configuration example of a fluid analysis apparatus according to the second embodiment of the present invention.

  The fluid analysis device 3 is built in, for example, the MRI device 4. However, the fluid analysis apparatus 3 may be connected to an MRI apparatus 4 or an image processing apparatus or an image server that stores image data collected in the MRI apparatus 4 via a network.

  The MRI apparatus 4 includes an imaging unit 4a. The imaging unit 4a has a function of collecting non-contrast fluid MR image data or contrast fluid MR image data for a plurality of frames by executing non-contrast imaging or contrast imaging of the fluid of the subject. Examples of the fluid include blood flow in blood vessels and cerebrospinal fluid (CSF). Imaging of blood vessels using an MRI apparatus is called magnetic resonance angiography (MRA).

Specifically, the imaging unit 4a includes a static magnetic field magnet, a gradient magnetic field coil, a gradient magnetic field power supply, a transmission radio frequency (RF) coil, a reception RF coil, a transmitter, a receiver, a sequencer, a computer, and the like. Has an element. In the imaging unit 4a, in accordance with imaging conditions such as a pulse sequence set in the computer, a gradient magnetic field is applied from the gradient magnetic field coil to the subject set under the static magnetic field in the static magnetic field magnet, and RF is transmitted from the transmission coil. Each pulse is applied. Then, a nuclear magnetic resonance (NMR) signal generated from the subject is received by the receiving coil, and an image reconstruction process is performed on the digitized NMR signal in the receiver by the computer, so that the MR image is obtained. Data is generated. Further, the gradient magnetic field pulse and the RF pulse are applied to the gradient magnetic field coil and the transmission RF coil from the gradient magnetic field power source and the transmitter, respectively, according to the pulse sequence output from the computer through the sequencer. Thereby, a gradient magnetic field pulse and an RF pulse having a desired waveform are applied to the subject from the gradient magnetic field coil and the transmission RF coil.
FIG. 7 is a diagram showing a first example of imaging conditions for non-contrast fluid imaging executed in the imaging unit 4a of the MRI apparatus 4 shown in FIG.

  In FIG. 7, the horizontal axis indicates time. FIG. 7 shows an example of a pulse sequence for rendering a non-contrast CSF with no periodicity. As shown in FIG. 7, a labeling pulse (also referred to as a tagging pulse) is applied to the CSF or the region of interest, followed by multiple acquisitions of imaging data (DATA ACQUISITION 1, DATA ACQUISITION 2, DATA ACQUISITION 3,...), It is possible to collect time-series non-contrast fluid image data for a plurality of frames corresponding to different times t1, t2, t3,.

  Labeling pulses include t-SLIP (Time Spatial Labeling Inversion Pulse), saturation (SAT) pulse (also called Rest grid pulse), SPAMM (spatial modulation of magnetization) pulse and DANTE pulse. Are known. The same type of labeling pulse may be applied a plurality of times, or a plurality of different kinds of labeling pulses may be applied in combination.

  For example, when a plurality of region selective 90 ° SAT pulses are combined, a region saturated with a radial or stripe pattern can be formed in the selected slab region. A region saturated with a desired pattern such as a stripe pattern, a grid pattern (lattice pattern), or a radial pattern can also be formed by applying a SPAMM pulse or a dante pulse. Therefore, it is possible to collect dynamic CSF image data for a plurality of frames in which the pattern moves with the CSF flow.

  In addition, imaging using t-SLIP includes a flow-in method and a flow-out method. The flow-in method is a method in which an area-selective inversion recovery (IR) pulse that constitutes t-SLIP is applied to the region of interest, and unlabeled CSF flowing into the region of interest from outside the region of interest is depicted. is there. On the other hand, in the flow-out method, a region non-selective IR pulse constituting t-SLIP is applied and a region selective IR pulse is applied to the labeling region, and the labeled CSF flowing from the labeling region to the region of interest is depicted. It is.

  The labeling pulse may be applied in synchronization with a biological signal such as an ECG (electro cardiogram) signal, a respiration sensor synchronization signal, or a pulse wave synchronization (PPG) signal.

On the other hand, as a sequence for collecting imaging data, an arbitrary sequence such as an SSFP (steady state free precession) sequence or a FASE (Fast ASE: fast asymmetric spin echo or fast advanced spin echo) sequence can be used.
FIG. 8 is a diagram illustrating a second example of imaging conditions for non-contrast fluid imaging executed in the imaging unit 4a of the MRI apparatus 4 shown in FIG.

  In FIG. 8, the horizontal axis indicates time, and the vertical direction also indicates time. FIG. 8 shows an example of a pulse sequence for imaging a subject's blood flow or CSF fluid using non-contrast imaging using IR pulses. That is, while gradually changing the inversion time (TI: inversion time) from the IR pulse application time to the imaging data acquisition timing to TI1, TI2, TI3, ..., the IR pulse application and imaging data acquisition (DATA ACQUISITION 1, DATA ACQUISITION 2, DATA ACQUISITION 3, ...) are executed repeatedly.

  The IR pulse is not limited to a single pulse, and a plurality of different or similar IR pulses may be applied at sufficiently short intervals that can be regarded as simultaneous. The region-selective IR pulse is applied to the imaging region in the case of the flow-in method described above, and is applied to the upstream region of the imaging region in the case of the flow-out method. The imaging region or fluid is then labeled by tilting the magnetization of the desired object by a combination of region selective IR pulses, frequency selective IR pulses, non-region selective IR pulses and non-frequency selective IR pulses. Then, fluid image data in which the fluid flowing into the imaging region is selectively depicted can be generated by image reconstruction processing of the imaging data collected from the imaging region at the timing when the fluid moves by a distance corresponding to TI. .

  Thereby, non-contrast fluid image data for a plurality of frames corresponding to a plurality of different TIs can be collected. In other words, when multiple frames of image data corresponding to multiple TIs are arranged in the shortest TI order, a continuous image is drawn in which fluid gradually flows as TI increases, such as dynamic image data. Data can be generated.

  The application of the IR pulse and the collection of imaging data may be synchronized with a biological signal such as an ECG signal in order to make the fluid velocity at the data collection timing closer to a constant value. Synchronous imaging using a biological signal is suitable for imaging fluid with periodicity such as blood flow.

As a pulse sequence for collecting imaging data, a pulse sequence such as an SE (spin echo) sequence, an FSE (fast spin echo) sequence, a FASE sequence, or an EPI (echo planar imaging) sequence is used.
FIG. 9 is a diagram showing a third example of imaging conditions for non-contrast fluid imaging executed in the imaging unit 4a of the MRI apparatus 4 shown in FIG.

  The application of the region selective IR pulse and the collection of the imaging data from the imaging region shown by the solid line including the blood flow as shown in FIG. 9 are repeated, and the application region of the region selective IR pulse is shown by the dotted line while keeping TI constant. In this way, continuous image data can be generated in which the fluid gradually flows like dynamic image data even by gradually changing to TAGGED REGION 1, TAGGED REGION 2, TAGGED REGION 3, ... be able to. In other words, the blood flow that existed in the area where the area-selective IR pulse was applied is labeled, and imaging data is collected at the timing when it flows out from the area-selective IR pulse application area downstream by a distance corresponding to TI. . Therefore, it is possible to generate image data in which the labeled blood flow is selectively depicted.

  Therefore, by collecting the imaging data repeatedly by changing the application area of the region-selective IR pulse, it corresponds to a plurality of different labeling regions, and the labeled blood flow is for a plurality of frames existing at different positions in the imaging region. Non-contrast fluid image data can be collected.

Furthermore, by combining the image data in which the blood flow is depicted upstream of the blood vessel with the image data corresponding to each labeling region through a synthesis process such as a maximum intensity projection (MIP) process, the blood flow gradually It is possible to generate continuous image data in which a state of flowing in the image is depicted.
In FIGS. 8 and 9, a SAT pulse may be used instead of the IR pulse.

  Up to this point, an example of imaging conditions for non-contrast imaging has been shown, but contrast image data can also be collected by a known contrast imaging method as in the first embodiment. In this case, contrast image data for a plurality of time-series frames corresponding to a plurality of time phases is collected.

  On the other hand, the fluid analysis device 3 has a function of acquiring information representing fluid dynamics based on non-contrast fluid MR image data or contrast fluid MR image data collected in the imaging unit 4a. For example, the fluid analysis device 3 causes the computer to read the fluid analysis program, thereby causing the computer to read the image data acquisition unit 3a, the tissue pixel change curve calculation unit 3b, the fluid pixel change curve calculation unit 3c, and the system matrix creation unit. It is made to function as 3d and the deconvolution part 3e. However, a part or all of the fluid analyzing apparatus 3 may be configured using a circuit.

  The image data acquisition unit 3a has a function of acquiring continuous fluid image data for a plurality of frames collected by the imaging unit 4a and supplying the fluid image data to the tissue pixel change curve calculation unit 3b.

  The tissue pixel change curve calculation unit 3b is configured such that, in continuous fluid image data for a plurality of frames, a change in pixel value between frames in a plurality of pixels corresponding to the tissue, a labeled blood flow that passes through the tissue, and the like. A function for calculating a change in the amount according to the concentration of the fluid and calculating each as a pixel density change curve of the tissue, and the calculated pixel density change curves of the plurality of tissues to the fluid pixel change curve calculating unit 3c and the system matrix creating unit 3d Has a function to output.

  In contrast, in the non-contrast imaging method, the contrast agent is not actually injected. For this reason, in the case of the non-contrast imaging method, a logical concentration corresponding to the ratio of the labeled fluid component to the total fluid amount is required, not the concentration related to the amount of fluid components such as blood flow.

For example, frames of time (t1, t2, t3, ...), TI (TI1, TI2, TI3, ...) or labeling areas (TAGGED REGION 1, TAGGED REGION 2, TAGGED REGION 3, ...) the pixel values corresponding to the i-th value of the parameter to be changed between the S _i, the initial value (t1, TI1, TAGGED REGION 1 ) pixel value corresponding to the value of the parameter as a reference, such as the S _0. Then, the change in the amount proportional to the pixel density in the tissue can be calculated as the pixel density change curve by the equation (10) shown in the first embodiment. Alternatively, a time intensity curve (TIC) representing a change in the pixel value S of the tissue corresponding to a parameter that can be changed between frames may be used as a pixel density change curve.
The amount and pixel value proportional to the pixel density are amounts corresponding to the fluid flow. Therefore, the pixel density change curve in the tissue includes fluid dynamic information in the tissue.

  The fluid pixel change curve calculation unit 3c calculates a pixel concentration change curve of a fluid circulating in a fluid region such as an artery based on the pixel concentration change curves in a plurality of pixels of the tissue, and the calculated pixel concentration of the fluid It has a function of outputting a change curve to the system matrix creation unit 3d. The fluid pixel density change curve is a curve representing a change between frames in an amount corresponding to the density of the pixel value corresponding to the fluid circulating in the fluid region.

  The pixel density change curve of the fluid can be calculated by the same method as in the first embodiment. Therefore, as in the first embodiment, the pixel density curve of the tissue is obtained for each of a plurality of pixels, whereas the pixel density curve of the fluid is a value representing the amount corresponding to the density of the pixel value of the fluid. Is obtained as a single curve. Therefore, when a region of interest (ROI) including a plurality of pixels is set as a fluid part, a value representing the amount corresponding to the density of the pixel value is calculated as in the first embodiment. The As representative values, there are a method using an average value between pixels and a method using a value calculated by curve fitting using a function of a constant multiple of a gamma probability density function having a fitting parameter as shown in Expression (11). . The pixels corresponding to the fluid region can also be specified automatically or by designation from the user as in the first embodiment.

  Based on the pixel density change curve corresponding to each pixel of the tissue, the system matrix creation unit 3d performs time (t1, t2, t3,...), TI (TI1, TI2, TI3,...) Or labeling region. A function of creating a one-dimensional column vector C having elements corresponding to the density of the pixel value of the tissue for each parameter value such as (TAGGED REGION 1, TAGGED REGION 2, TAGGED REGION 3, ...); And a function of creating a one-dimensional column vector whose element is an amount corresponding to the density of the pixel value of the fluid for each parameter value based on the pixel density change curve of the fluid. In addition, the system matrix creating unit 3d creates a C = AX system matrix using a column vector C corresponding to the tissue, a coefficient matrix A in which the column vectors corresponding to the fluid are two-dimensionally arranged, and an unknown quantity X. And a function of outputting the created system matrix to the deconvolution unit 3e.

The column vector C corresponding to the tissue, the column vector corresponding to the fluid, the coefficient matrix A, and the system matrix C = AX can be created by a method similar to the method in the first embodiment. In other words, the description is as follows.
FIG. 10 is a diagram for explaining a method of creating an element to be analyzed in the system matrix creation unit 3d shown in FIG.

  In FIG. 10, the horizontal axis represents time (t1, t2, t3,...), TI (TI1, TI2, TI3,...) Or labeling region (TAGGED REGION 1, TAGGED REGION 2, TAGGED REGION 3,. .) And the like, and the vertical axis indicates the density of the pixel value.

As shown in FIG. 10, when pixel density curves C1 _i , C2 _i , C3 _i ,... And fluid pixel density curves a _i corresponding to each pixel of the tissue are obtained, a plurality of tissue pixels. A first point before at least rising of the density change curves C1 _i , C2 _i , C3 _i ,... And a second point before at least the rising of the pixel density change curve a _{i of the} fluid are determined. For example, the first point is the time (t1, t2, t3,...), TI (TI1, TI2, TI3,...) Or the labeling region (TAGGED REGION 1, TAGGED REGION 2, TAGGED REGION 3,. ..) and the like are changed to initial values (t1, TI1, TAGGED REGION 1) of parameters changed in the imaging. In addition, the second point is determined to have a value later than the first point.

Furthermore, the analysis target element corresponding to the tissue pixel density curve C1 _i , C2 _i , C3 _i , ... and the analysis target element corresponding to the fluid pixel density curve a _i are the same number. , The analysis target ranges RANGE1, RANGE2 having the same width starting from the first point and the second point are the pixel density curves C1 _i , C2 _i , C3 _i ,. Set for each _i . It should be noted that some or all of the elements in front of the parameter values included only in the analysis target range RANGE1 of the tissue pixel density curves C1 _i , C2 _i , C3 _i ,. On the other hand, some or all of the elements behind the parameter values included only in the analysis target range RANGE2 of the fluid pixel density curve a _i may be zero.

As a result, the amount or zero according to the density of the tissue pixel value at the first point is the first element to be analyzed corresponding to the tissue pixel density curve C1 _i , C2 _i , C3 _i ,. The amount corresponding to the density of the pixel value of the fluid at the second point is the first element to be analyzed corresponding to the pixel density curve a _i of the fluid.

  Using the elements corresponding to the pixel density curve of the same number of tissues and the elements corresponding to the pixel density curve of the fluid thus created, a system matrix as shown in Expression (15) can be created.

  The deconvolution unit 3e has the same function as the function of the deconvolution unit 2e in the first embodiment. That is, the deconvolution unit 3e has a function of calculating the unknown number X by deconvolution integration of a coefficient matrix A in which a column vector C corresponding to the tissue and a column vector corresponding to the fluid are two-dimensionally arranged, and an impulse obtained from the unknown number X. Based on the response function, it has a function of calculating information on fluid dynamics such as an average fluid transit time MTT, local cerebral blood flow CBF, local cerebral blood flow CBV, local CSF flow rate, and local CSF flow.

Moreover, the deconvolution part 3e is comprised so that it may correct | amend the information regarding fluid dynamics as shown in Formula (19) similarly to the deconvolution part 2e in 1st Embodiment as needed.
Next, the operation and action of the fluid analysis device 3 will be described.

FIG. 11 is a flowchart showing the flow of data processing in the fluid analyzing apparatus 3 shown in FIG. Note that the same processing as the detailed processing shown in FIG. 4 in the first embodiment is omitted.
First, the imaging unit 4a of the MRI apparatus 4 performs non-contrast imaging or contrast imaging in advance and collects MR fluid image data for a plurality of frames.

  In step S21, the image data acquisition unit 3a of the fluid analysis device 3 acquires fluid image data for a plurality of frames from the imaging unit 4a and supplies the fluid image data to the tissue pixel change curve calculation unit 3b.

  Next, in step S22, the tissue pixel change curve calculation unit 3b calculates the change in the amount according to the density of the pixel value between frames in a plurality of pixels corresponding to the tissue of the fluid image data for a plurality of frames. The pixel density change curve is calculated.

  Next, in step S23, the fluid pixel change curve calculation unit 3c calculates a fluid pixel concentration change curve based on the pixel concentration change curves in a plurality of pixels of the tissue.

  Next, in step S24, the system matrix creation unit 3d corresponds the column vector corresponding to the pixel density of the tissue to C and the pixel density of the fluid based on the pixel density change curve of the tissue and the pixel density change curve of the fluid. A system matrix of C = AX is created in which a coefficient matrix in which column vectors are two-dimensionally arranged is A and an unknown is X. At this time, the first element of the column vector C corresponding to the tissue is a value corresponding to the first parameter value on the horizontal axis of the pixel density change curve of the tissue. On the other hand, the first element of the column vector corresponding to the fluid is a value corresponding to the second parameter value behind the first parameter value on the horizontal axis of the pixel density change curve of the fluid. Further, each column vector is created so that the number of elements of the column vector C corresponding to the tissue is the same as the number of elements of the column vector corresponding to the fluid. Therefore, zero may be used instead of pixel density as an element constituting the column vector.

  Next, in step S25, the deconvolution unit 3e calculates an unknown number X of the system matrix, and calculates information related to fluid dynamics based on an impulse response function obtained from the unknown number X. This fluid dynamic information may be corrected by equation (19).

  Thus, the fluid analysis device 3 is not limited to time-series contrast image data corresponding to a plurality of time phases, but based on non-contrast image data for a plurality of continuous frames collected by changing imaging parameters. Information on fluid dynamics such as blood flow and CSF can be obtained by processing similar to that for contrast image data.

  In the second embodiment, the example in which the fluid analysis apparatus 3 acquires information on fluid dynamics based on image data collected by the MRI apparatus has been described. However, an X-ray CT apparatus, an ultrasonic diagnostic apparatus, or a PET apparatus is described. It is also possible to configure the fluid analysis device 3 to acquire information on fluid dynamics based on contrast image data or non-contrast image data collected by a medical image capturing device such as the above.

(Other embodiments)
Although specific embodiments have been described above, the described embodiments are merely examples, and do not limit the scope of the invention. The novel methods and apparatus described herein can be implemented in a variety of other ways. Various omissions, substitutions, and changes can be made in the method and apparatus described herein without departing from the spirit of the invention. The appended claims and their equivalents include such various forms and modifications as are encompassed by the scope and spirit of the invention.

DESCRIPTION OF SYMBOLS 1 Blood flow dynamic analysis apparatus 2 Local blood flow analysis application 2a Local blood flow analysis application 2b Contrast density curve calculation part 2c Arterial time density curve calculation part 2d System matrix creation part 2e Deconvolution part 3 Fluid analysis apparatus 3a Image data acquisition part 3b Tissue pixel change curve calculation unit 3c Fluid pixel change curve calculation unit 3d System matrix creation unit 3e Deconvolution unit 4 MRI apparatus 4a Imaging unit

Claims (13)

An amount corresponding to the concentration of the contrast agent in the tissue by measuring the change in the pixel value of each pixel between the plurality of time phases based on the contrast image data of the subject in the plurality of time phases obtained by the medical imaging apparatus A tissue time concentration curve calculating means for calculating a time concentration curve representing a time change of
Arterial time concentration curve calculating means for setting an arterial region from the time concentration curve corresponding to the tissue and calculating a time concentration curve corresponding to the artery;
The time density curve corresponding to the tissue is discretized with reference to a first time, and a plurality of elements including values corresponding to the concentration of the contrast agent in the tissue at different times are arranged in one dimension. First constructing means constituting one set;
The time concentration curve corresponding to the artery is discretized with reference to a second time behind the first time, and includes a plurality of values corresponding to the concentration of the contrast agent in the artery at different times A second constructing means constituting a second set in which the elements are arranged two-dimensionally;
Deconvolution integration means for performing deconvolution of the first set and the second set and calculating a transfer function of the tissue;
A blood flow dynamics analysis device comprising: blood flow information calculation means for calculating information related to blood flow dynamics based on the transfer function. The first time is a collection start time of the image data or a predetermined time between the collection start time and the second time, and the second time is the time concentration curve corresponding to the artery. The blood flow dynamic analysis device according to claim 1, which is a time at which the blood pressure starts to rise. The first set includes H first sets in a time axis direction with respect to a set of values corresponding to the concentration of the contrast agent in the tissue at L different times arranged in one dimension. The second set is a time axis relative to a set of values corresponding to the contrast agent concentration in the artery at L different times arranged in two dimensions. The blood flow dynamics analysis device according to claim 1 or 2, which is a set in which H second elements are added to the rear of the direction. The first element is at least one of 0 and a value corresponding to the concentration of the contrast agent in the tissue, and the second element is at least one of a value corresponding to 0 and the concentration of the contrast agent in the artery. The blood flow dynamics analysis device according to claim 3 which is one side. The blood flow information calculation means searches for the maximum value of the transfer function, calculates a peak range from the searched maximum value of the transfer function, and based on the maximum value of the transfer function and the peak range, the blood flow The blood flow dynamics analysis device according to any one of claims 1 to 4, which calculates information relating to dynamics. Computer
An amount corresponding to the concentration of the contrast agent in the tissue by measuring the change in the pixel value of each pixel between the plurality of time phases based on the contrast image data of the subject in the plurality of time phases obtained by the medical imaging apparatus A tissue time concentration curve calculating means for calculating a time concentration curve representing a time change of
An arterial time concentration curve calculating means for setting an arterial region from the time concentration curve corresponding to the tissue and calculating a time concentration curve corresponding to the artery;
The time density curve corresponding to the tissue is discretized with reference to a first time, and a plurality of elements including values corresponding to the concentration of the contrast agent in the tissue at different times are arranged in one dimension. First constituting means constituting one set;
The time concentration curve corresponding to the artery is discretized with reference to a second time behind the first time, and includes a plurality of values corresponding to the concentration of the contrast agent in the artery at different times A second constituting means constituting a second set in which the elements are arranged two-dimensionally;
Deconvolution integration means for performing deconvolution integration of the first set and the second set and calculating a transfer function of the tissue, and blood flow information for calculating information related to blood flow dynamics based on the transfer function Calculation means,
Blood flow dynamics analysis program to function as. Based on the image data of the plurality of frames of the subject, the first change between the plurality of frames is calculated in an amount corresponding to the density of the pixel value in the pixel corresponding to the tissue, while the pixel value corresponding to the fluid is calculated. Density change calculating means for calculating a second change between the plurality of frames in an amount corresponding to the density;
Based on the first change, a first plurality of elements to be analyzed are created with the element corresponding to the first point of the first change as the first element, while the second change The second plurality of elements to be analyzed with the element corresponding to the second point after the first point of the second change as the first element based on Fluid information calculation means for calculating information on fluid dynamics with deconvolution integration processing based on the first plurality of elements and the second plurality of elements;
A fluid analysis apparatus comprising: 8. The density change calculating unit is configured to calculate a time change of an amount corresponding to a concentration of a contrast agent as the first and second changes based on contrast image data of a plurality of frames. The fluid analysis apparatus described. The fluid analysis apparatus according to claim 7, wherein the concentration change calculation unit is configured to calculate the first and second changes based on non-contrast image data of a plurality of frames. The density change calculating means is configured to change the time from the common labeling pulse to the acquisition timing of the imaging data by the magnetic resonance imaging apparatus based on the non-contrast image data of a plurality of frames dynamically acquired. The fluid analyzing apparatus according to claim 9, wherein the fluid analyzing apparatus is configured to calculate a change in the pressure. The concentration change calculation means repeats the application of the inversion recovery pulse or saturation pulse and the collection of imaging data by the magnetic resonance imaging apparatus, and the time from the application timing of the inversion recovery pulse or saturation pulse to the imaging data collection timing The fluid analysis device according to claim 9, configured to calculate the first and second changes based on non-contrast image data of a plurality of frames collected by changing the frequency. The density change calculation means repeats the application of the inversion recovery pulse or saturation pulse and the collection of imaging data by the magnetic resonance imaging apparatus, and changes the application area of the inversion recovery pulse or saturation pulse, and acquires a plurality of frames collected. The fluid analysis apparatus according to claim 9, configured to calculate the first and second changes based on non-contrast image data. Computer
Based on the image data of the plurality of frames of the subject, the first change between the plurality of frames is calculated in an amount corresponding to the density of the pixel value in the pixel corresponding to the tissue, while the pixel value corresponding to the fluid is calculated. A density change calculating means for calculating a second change between the plurality of frames in an amount corresponding to the density; and an element corresponding to the first point of the first change based on the first change While the first plurality of elements to be analyzed as elements are created, the second change corresponds to a second point after the first point of the second change based on the second change A second plurality of elements to be analyzed, the number of which is the first element, is created, and a fluid with deconvolution integration processing is performed based on the first plurality of elements and the second plurality of elements. Information calculation means for calculating information on the dynamics of the fluid ,
Fluid analysis program to function as

Images