WO2021095868A1 - 評価装置、評価方法およびプログラム - Google Patents

評価装置、評価方法およびプログラム Download PDF

Info

Publication number
WO2021095868A1
WO2021095868A1 PCT/JP2020/042497 JP2020042497W WO2021095868A1 WO 2021095868 A1 WO2021095868 A1 WO 2021095868A1 JP 2020042497 W JP2020042497 W JP 2020042497W WO 2021095868 A1 WO2021095868 A1 WO 2021095868A1
Authority
WO
WIPO (PCT)
Prior art keywords
signal
time
value
polarization
measurement signal
Prior art date
Application number
PCT/JP2020/042497
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 US17/776,191 priority Critical patent/US20220390357A1/en
Priority to JP2021556189A priority patent/JPWO2021095868A1/ja
Publication of WO2021095868A1 publication Critical patent/WO2021095868A1/ja

Links

Images

Classifications

    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/21Polarisation-affecting properties
    • G01N21/23Bi-refringence
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/48Biological material, e.g. blood, urine; Haemocytometers
    • G01N33/483Physical analysis of biological material
    • G01N33/4833Physical analysis of biological material of solid biological material, e.g. tissue samples, cell cultures
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N21/47Scattering, i.e. diffuse reflection
    • G01N21/4795Scattering, i.e. diffuse reflection spatially resolved investigating of object in scattering medium
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/48Biological material, e.g. blood, urine; Haemocytometers
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N33/00Investigating or analysing materials by specific methods not covered by groups G01N1/00 - G01N31/00
    • G01N33/48Biological material, e.g. blood, urine; Haemocytometers
    • G01N33/483Physical analysis of biological material
    • GPHYSICS
    • G01MEASURING; TESTING
    • G01NINVESTIGATING OR ANALYSING MATERIALS BY DETERMINING THEIR CHEMICAL OR PHYSICAL PROPERTIES
    • G01N21/00Investigating or analysing materials by the use of optical means, i.e. using sub-millimetre waves, infrared, visible or ultraviolet light
    • G01N21/17Systems in which incident light is modified in accordance with the properties of the material investigated
    • G01N2021/178Methods for obtaining spatial resolution of the property being measured
    • G01N2021/1785Three dimensional
    • G01N2021/1787Tomographic, i.e. computerised reconstruction from projective measurements

Definitions

  • the present invention relates to evaluation devices, methods and programs.
  • the present invention relates to a technique for visualizing and quantitatively evaluating the state of a sample such as a living tissue by processing a measurement signal measured by an optical coherence tomography (OCT).
  • OCT optical coherence tomography
  • Non-Patent Document 1 a signal analysis method called “dynamic OCT” has been proposed (Non-Patent Document 1).
  • this method is poorly quantitative, and it is difficult to correctly evaluate the degree of activity of a living body by this method.
  • this method is suitable for a special type of OCT called Full-field OCT (FF-OCT), and it is difficult to implement it in a widely used scanning OCT.
  • FF-OCT Full-field OCT
  • the present invention has been made in view of the above circumstances, and it is an object of the present invention to provide an evaluation method capable of quantitatively evaluating the dynamic characteristics of a sample, for example, the dynamics of living tissue and intracellular activity. Make one.
  • the present invention employs the following means.
  • One aspect of the present invention is to acquire an optical coherence tomography (OCT) signal indicating the state of a living tissue as a sample, and obtain a signal value based on the OCT signal for each observation point in the sample. It is an evaluation device having a measurement unit to be acquired and an evaluation unit for calculating a time fluctuation characteristic value indicating a time fluctuation characteristic of the signal value within a predetermined period.
  • OCT optical coherence tomography
  • the evaluation unit may calculate the variance of the signal value as the time fluctuation characteristic value.
  • the evaluation unit sets the sum of squares of the deviations between the signal strength of the OCT signal and the average value of the signal strength at each frame time within the predetermined period as the frame in the predetermined period.
  • the variance may be calculated for each observation point by dividing by a number.
  • the evaluation unit calculates a correlation coefficient between the signal value and a signal value obtained by time-shifting the signal value by a time shift amount ⁇ for each time shift amount ⁇ .
  • the decay rate of the correlation coefficient with the increase of the time shift amount may be calculated as the time fluctuation characteristic value.
  • the evaluation unit calculates the sum of squares of the deviations between the signal strength of the OCT signal and the average value of the signal strength at each frame time within the predetermined period as a variance.
  • the sum of the products of the time-shifted signal strength and the deviation from the average value is calculated as the covariance, and the covariance is divided by the variance to be calculated as the correlation coefficient for each time shift amount ⁇ . Even if regression analysis is performed using a predetermined attenuation function using the correlation coefficient for each time shift amount ⁇ , and the parameter of the attenuation function that approximates the correlation coefficient is calculated as the attenuation rate for each observation point. Good.
  • the evaluation unit may calculate the attenuation rate by using the correlation coefficient calculated with the time shift amount ⁇ as non-zero.
  • the first incident component interferes with the component reflected or scattered from the sample with respect to the first incident component incident on the sample in the first polarized state.
  • the first measurement signal having the first polarization state
  • the second measurement signal having the second polarization state with respect to the first interference component
  • the second incident component incident on the sample in the second polarization state It has a third measurement signal having a first polarization state among the second interference components in which the second incident component interferes with a component reflected or scattered from the sample, and a second polarization state with respect to the second interference component.
  • the polarization characteristic value based on the polarization characteristic at the observation point in the sample may be determined, and the evaluation unit may determine the time fluctuation characteristic value indicating the time variation characteristic of the polarization characteristic value. ..
  • the measuring unit uses a Jones matrix for each observation point based on the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal.
  • the cumulative Jones matrix at the observation point is determined from the Jones matrix at the observation point in the sample and the Jones matrix at the surface of the sample, and the polarization characteristic value is the cumulative phase difference between the eigenvalues of the cumulative Jones matrix.
  • a phase lag index value may be set.
  • the measuring unit uses a Jones matrix for each observation point based on the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal.
  • the local Jones matrix between the Jones matrix at the first observation point in the sample and the Jones matrix at the second observation point in the sample to the first observation point and the second observation point is determined, and the local The polarization characteristic value may be determined based on the local phase delay, which is the phase difference between the eigenvalues of the Jones matrix.
  • the measuring unit divides the local phase delay by the wave number of the incident light incident on the sample and the thickness of the first observation point and the second observation point.
  • the birefringence may be determined.
  • the evaluation unit may calculate the time variation characteristic value based on the dispersion or standard deviation of the polarization characteristic value.
  • the evaluation unit may calculate the time variation characteristic value based on the logarithmic variance or standard deviation of the polarization characteristic value.
  • the evaluation unit may calculate the dynamic contrast by dividing the standard deviation of the polarization characteristic value by the average value of the birefringence.
  • the measuring unit uses the first measurement signal, the first Jones vector based on the second measurement signal, the third measurement signal, and the fourth measurement signal as the polarization characteristic value.
  • the second Jones vector based on the measurement signal is converted into the first Stokes vector and the second Stokes vector, respectively, and the evaluation unit uses the time average value of the first Stokes vector and the second Stokes vector as the time fluctuation characteristic values.
  • the time polarization uniformity may be determined based on the time average value of the Stokes vector.
  • the measuring unit has the time average value of the corrected first Stokes vector obtained by subtracting the noise component from the first Stokes vector, and the second Stokes vector.
  • the time polarization uniformity may be determined based on the time average value of the corrected second Stokes vector obtained by subtracting the noise component from.
  • the measuring unit uses the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal as the polarization characteristic values.
  • a Jones matrix may be determined for each observation point, and the evaluation unit may calculate the von Neumann entropy of the Jones matrix as the time fluctuation characteristic value.
  • the evaluation unit uses the first measurement signal and the first Jones vector based on the second measurement signal, and the third measurement signal and the second measurement signal based on the fourth measurement signal.
  • the entropy of the noise component is calculated from the time polarization uniformity of the converted first Stokes vector and the time polarization uniformity of the second Stokes vector, respectively, and the von Neumann entropy is based on the entropy of the noise component. It may be corrected.
  • the first polarized state is horizontally polarized light
  • the second polarized state is vertically polarized light
  • the first measurement signal is a first horizontally polarized light spectrum interference signal.
  • the second measurement signal is a second horizontally polarized spectrum interference signal
  • the third measurement signal is a first vertically polarized spectrum interference signal
  • the fourth measurement signal is a second vertically polarized spectrum interference signal.
  • the evaluation unit may calculate the time fluctuation characteristic value for each observation period interval longer than the predetermined period.
  • Another aspect of the present invention may include an output processing unit that determines an evaluation value indicating the active state of the sample based on the time fluctuation characteristic value.
  • image data having an output value for the time fluctuation characteristic value for each observation point as a signal value by using a function that gives an output value that changes monotonically with respect to a change in an input value.
  • An image processing unit may be provided.
  • Another aspect of the present invention is a method in an evaluation device, in which an optical coherence tomography (OCT) signal indicating the state of a living tissue as a sample is acquired and a signal value based on the OCT signal is obtained.
  • OCT optical coherence tomography
  • an optical coherence tomography (OCT) signal indicating the state of a living tissue as a sample is acquired by a computer of an evaluation device, and a signal value based on the OCT signal is obtained.
  • OCT optical coherence tomography
  • TPU time change of TPU. It is a figure which shows an example of the time change of the dynamic contrast of a birefringence index. It is a figure which shows an example of the time change of the logarithmic intensity dispersion.
  • FIG. 1 is a diagram illustrating a method for evaluating a living tissue according to the first embodiment of the present invention.
  • the method for evaluating the living tissue of the present embodiment is that the OCT system 1 includes an imaging means (imaging unit) 10 that performs OCT (optical coherence tomography) imaging of the sample Sm a plurality of times within a predetermined time, and each time.
  • the living tissue evaluation device 20 provided with the above, the minute fluctuations of the living tissue are visualized and quantitatively evaluated according to the following procedure (see FIG. 8).
  • OCT optical coherence tomography
  • measurement is performed (several hundred frames from the frame), preferably 10 frames or more, more preferably 15 frames or more.
  • the predetermined period (observation period) related to one OCT imaging is shorter than the observation period interval until the adjacent observation period, and the dynamic characteristics of the sample appearing in the OCT signal obtained by the OCT imaging at that time are evaluated. Any period may be sufficient as long as the required accuracy can be ensured.
  • the observation period is set between the frame interval so that the frequency component of the dynamic characteristics of the sample is included in the frequency band from the lowest frequency corresponding to the observation period to the highest frequency corresponding to the frame interval.
  • the observation period interval may be a period during which the accuracy required for the evaluation of the global change tendency of the activity state of the biological tissue as a sample can be ensured. For example, a series of processes (apoptosis) from the mother's body to the sample living tissue being separated from the mother's body to cell death, or a series of processes required to cell death after a cause occurs in the sample's body tissue in the mother's body. A time sufficiently shorter than the required time such as (necrosis) may be set as the observation period interval. From the viewpoint of knowing the characteristics of the entire sample, it is preferable to perform the same OCT imaging on other parts.
  • the measuring means 22 and the evaluation means 24 quantitatively evaluate the activity of the living tissue based on the time change of the OCT signal intensity obtained from each OCT image.
  • Specific methods for quantitative evaluation include calculation of speckle variance (SV: Speckle Variance), attenuation speed of OCT correlation (OCDS: OCT Correlation Decay Speed), and evaluation of the calculation result.
  • the speckle variance indicates the dispersion (fluctuation) of the OCT signal strength in a short time, and can be calculated by using the following equation (1).
  • x, z is the surface of each biological tissue indicates the position in the depth direction from the surface
  • I (x, z, t i) is displayed on a linear scale or a logarithmic scale, etc.
  • the OCT signal strength at each position and each time is shown
  • ⁇ I> shows the average OCT signal strength
  • N shows the number of frames of the OCT signal within a predetermined period.
  • x and z correspond to the positions of the observation points corresponding to the individual pixels forming the OCT image.
  • the measuring means 22 first performs the OCT signal strength I (x, z, t) for each observation point (x, z) at the time t 1 , t 2 , ... T N of each frame. 1 ), I (x, z, t 2 ), ... I (x, z, t N ) are measured. Subsequently, the evaluation means 24 calculates the average ⁇ I> of the OCT signal intensities for N times. Subsequently, the evaluation means 24 squares the difference between the OCT signal strength of each time and the average ⁇ I>.
  • the evaluation means 24 is [I (x, z, t 1 )- ⁇ I>] 2 , [I (x, z, t 2 )- ⁇ I>] 2 , ... [I (x, z). , T N )- ⁇ I>] 2 is calculated. Finally, the evaluation means 24 can obtain the speckle variance ⁇ 2 (x, z) by dividing the sum of these by N.
  • the decay rate of the OCT correlation indicates the rate at which the correlation coefficient between adjacent times (t, t + ⁇ ) decreases as the time shift amount ⁇ increases.
  • ⁇ (x, z, ⁇ ) can be calculated.
  • x and z indicate the positions in the surface of the living tissue in the depth direction from the surface, respectively, and ⁇ cov 2 (x, z, ⁇ ) and ⁇ I 2 (x, z) are.
  • the co-dispersion and speckle variance (dispersion) of the OCT signal strengths I (x, z, t) and I (x, z, t + ⁇ ) are shown, respectively.
  • the OCT signal strength I (x, z, t + ⁇ ) is a signal value obtained by time-shifting the OCT signal strength I (x, z, t) with a time shift amount ⁇ . That is, the correlation coefficient ⁇ (x, z, ⁇ ) corresponds to the auto-correlation function of the signal value I (x, z, t) of the OCT signal at the observation point (x, z).
  • the evaluation means 24 sets the covariance ⁇ cov 2 (x, z, ⁇ ) between adjacent times (t, t + ⁇ ) as [I (x, z, t) ⁇ ⁇ I>] ⁇ . [I (x, z, t + ⁇ )- ⁇ I>] is calculated. Subsequently, the evaluation means 24 divides the covariance ⁇ cov 2 (x, z, ⁇ ) by the speckle variance ⁇ 2 (x, z) calculated by the above equation (1) to obtain a correlation coefficient ⁇ ( x, z, ⁇ ) can be obtained.
  • the evaluation means 24 applies the obtained correlation coefficient ⁇ (x, z, ⁇ ) to a predetermined function of the time shift amount ⁇ by performing regression analysis (regression), and the function value of that function.
  • a predetermined function may be a function in which the function value given to the time shift amount is attenuated as the time shift amount increases, for example, an exponential function. When an exponential function is used, its parameter bottom is obtained as an index of damping rate.
  • Linear analysis may be used as the method of regression analysis, but the method is not limited to this, and non-linear analysis may be used.
  • the method for evaluating a living tissue is to acquire an image of a specific part of the living tissue many times in a short time and calculate the speckle variance, the decay rate of the OCT correlation, and the like. Quantitatively evaluates the activity of living tissue. As a result, minute changes (fluctuations) that have been overlooked in the past can be detected, and quantitative evaluation of the dynamics of living tissues can be realized.
  • Example 1 100 high-speed OCT images were performed every 2 hours at 13 ms intervals on the same tumor aggregate sample as in Comparative Example 1.
  • the observation period interval was set to 2 hours in order to capture the tendency of changes in the active state during the period from the removal of the sample tumor aggregate tissue from the living body to death (typically about 1 to 3 days). Because it is enough.
  • the reason why one observation period was set to 100 times at 13 ms intervals, that is, 1.3 seconds, is that this observation period is sufficiently shorter than the observation period interval, and the movement or intracellular activity of the cells constituting the tissue. This is because it is sufficient to grasp the time change of the optical characteristics according to (for example, several Hz to 20 Hz).
  • the speckle variance was calculated based on the time change of the OCT signal intensity obtained from each OCT image.
  • Example 2 The same high-speed OCT imaging as in Example 1 was performed on the tumor aggregate sample similar to Comparative Example 1. Subsequently, the decay rate of the OCT correlation was calculated based on the time change of the OCT signal intensity obtained from each OCT image.
  • Example 3 The same high-speed OCT imaging as in Example 1 was performed on the tumor aggregate sample similar to Comparative Example 1. Subsequently, the attenuation coefficient (AC: Attenuation Coefficient) of the OCT signal intensity in the depth direction was calculated based on the OCT signal intensity obtained from each OCT image.
  • AC Attenuation Coefficient
  • FIG. 2 is an image of a tumor aggregate sample obtained by applying the methods of Comparative Examples 1 and 1 and 2.
  • the upper, middle, and lower images correspond to Comparative Examples 1 and 1, respectively.
  • the numerical values of 0 hr, 8 hr, and 28 hr indicate the elapsed time from the time when the culture was started after being cut out from the human body, respectively.
  • the state of the tumor aggregate sample does not change with the passage of time.
  • the images of Examples 1 and 2 there is a dark part in the center and a bright part around the dark part in the initial stage, but it can be seen that the bright part becomes darker with the passage of time. Since it is considered that the dark part represents the dead state and the bright part represents the living state, it is possible to confirm the life-and-death state of the tumor aggregate sample from the changes in these images.
  • FIGS. 3A and 3B are diagrams showing changes over time in the decay rate of speckle variance and OCT correlation calculated for each observation period from the samples shown in the images of Examples 1 and 2 of FIG. 2, respectively.
  • FIG. 3C is a diagram showing the attenuation coefficient of the OCT signal intensity obtained by applying the method of Example 3.
  • the horizontal axes of FIGS. 3A, 3B, and 3C all show the elapsed time since the tumor aggregate sample was cut out, and the vertical axis shows the speckle variance and OCT correlation of the entire tumor aggregate sample, respectively.
  • the average value of the decay rate and the attenuation coefficient of the OCT signal strength is shown. In each figure, a decreasing tendency of each average value can be seen with the passage of time.
  • the first of the plurality of OCT images does not correctly reflect the attenuation tendency of the correlation coefficient with the subsequent samples, and causes a significant difference from the estimated value estimated from the attenuation tendency.
  • the correlation coefficient when the time shift amount ⁇ becomes zero is excluded, and the decay rate of the OCT correlation is calculated using the correlation coefficient with respect to the non-zero time shift amount ⁇ which is non-zero. Is preferable.
  • FIGS. 4A and 4B are diagrams showing the life and death states of living tissues to which the methods of Examples 1 and 2 are applied, respectively.
  • the horizontal axis of FIGS. 4A and 4B shows the elapsed time, and the vertical axis shows the content of living cells or dead cells.
  • the boundary between life and death was set to 3.0 in Example 1 (SV) and 5 ⁇ 10 -4 ms -1 in Example 2 (OCDS). In each graph, there is a tendency that the number of living cells decreases and the number of dead cells increases with the passage of time.
  • Example 4 The same liver sample as in Comparative Example 2 was subjected to 100 high-speed OCT imaging at 10 ms intervals every hour. Subsequently, the speckle variance was calculated based on the time change of the OCT signal intensity obtained from each OCT image.
  • Example 5 The same high-speed OCT imaging as in Example 4 was performed on the same liver sample as in Comparative Example 2. Subsequently, the decay rate of the OCT correlation was calculated based on the time change of the OCT signal intensity obtained from each OCT image.
  • Example 6 The same high-speed OCT imaging as in Example 4 was performed on the same liver sample as in Comparative Example 2. Subsequently, the attenuation coefficient (AC) of the OCT signal intensity in the depth direction was calculated based on the OCT signal intensity obtained from each OCT image.
  • AC attenuation coefficient
  • FIG. 5 is an image of a liver sample obtained by applying the methods of Comparative Examples 2, 4 and 5.
  • the upper, middle, and lower images correspond to Comparative Example 2, Examples 4, and 5, respectively.
  • Numerical values such as 0 hr, 8 hr, and 16 hr indicate the elapsed time from the time when they were cut out from the mother mouse, respectively.
  • the state of the liver sample does not change with the passage of time.
  • the images of Examples 4 and 5 there is a dark part on the lower side and a bright part on the upper side in the initial stage, but it can be seen that the bright part becomes darker with the passage of time. Since it is considered that the dark part represents the dead state and the bright part represents the living state, it is possible to confirm the alive / dead state of the liver sample from the changes in these images.
  • FIGS. 6A and 6B are diagrams showing the decay rate of speckle variance and OCT correlation calculated for each observation period from the samples shown in the images of Examples 4 and 5, respectively.
  • FIG. 6C is a diagram showing the attenuation coefficient of the OCT signal intensity obtained by applying the method of Example 6.
  • the horizontal axes of FIGS. 6A, 6B, and 6C all show the elapsed time since the liver sample was cut out, and the vertical axis shows the speckle variance and the decay rate of the OCT correlation in the entire liver sample, respectively.
  • the average value of the attenuation coefficient is shown. In each figure, with the passage of time, a two-step decreasing tendency with different slopes can be seen for each average value.
  • the OCT correlation is most suitable for quantitative evaluation because the slope of the decay rate is large.
  • the magnitude of the slope with respect to the speckle variance is the next largest, and is significantly larger than the slope with respect to the attenuation coefficient of the OCT signal strength proposed conventionally.
  • the first image of the plurality of OCT images does not correctly reflect the attenuation tendency of the correlation coefficient with the subsequent samples, and causes a significant difference from the estimated value estimated from the attenuation tendency.
  • FIGS. 7A and 7B are diagrams showing the life and death states of living tissues to which the methods of Examples 4 and 5 are applied, respectively.
  • the horizontal axis of FIGS. 7A and 7B shows the elapsed time, and the vertical axis shows the content of living cells or dead cells.
  • the boundary between life and death was set to 3.0 in Example 4 (SV) and 5 ⁇ 10 -4 ms -1 in Example 5 (OCDS).
  • SV speckle variance
  • OCDS OCT correlation
  • FIG. 9 is a configuration diagram showing an example of the OCT system 1 according to the present embodiment.
  • the OCT system 1 constitutes PS-OCT.
  • the PS-OCT includes an optical system for irradiating the sample Sm with incident light having a known polarization state and acquiring interference light in which the reflected light reflected from the sample Sm and the reference light interfere with each other.
  • the OCT system 1 includes a measurement signal processing device 200 that analyzes the change characteristic of the polarization state in the sample Sm from the polarization state of the interference light acquired by the optical system.
  • the measurement signal processing device 200 functions as an evaluation device that analyzes the state of the living tissue used as the sample Sm using the OCT signal.
  • the measurement signal processing device 200 generates an image that visualizes the analyzed change characteristics.
  • the object to be observed as the sample Sm is mainly a part of a living body such as a human or an animal. More specifically, it may be any of fundus, blood vessels, teeth, subcutaneous tissue and the like. As a result, the state inside the sample Sm can be measured or observed non-invasively. Therefore, it is expected to be applied to the diagnosis of in vivo tissues, for example, microvessels such as the fundus, lymphatic vessels, and myocardium.
  • the OCT system 1 illustrated in FIG. 9 is an observation system that applies a wavelength sweep type OCT (SS-OCT: Swept Source-OCT) for sweeping the wavelength of light generated by the light source 102 to obtain a spectral interference signal. Nasu.
  • the OCT system 1 branches the light emitted from the light source 102 to the probe arm (described later) and the reference arm 130 to enter the light.
  • the OCT system 1 separates the light branched to the probe arm into a horizontally polarized light component and a vertically polarized light component, and scans (B-scan) the light containing the polarized light components having different optical path lengths between them on the sample Sm to be measured.
  • the OCT system 1 acquires the interference light by interfering the reference light branched to the reference arm 130 with the reflected light which is a component obtained by reflection, scattering, or both of the sample Sm.
  • the direction in which light is irradiated to the sample Sm is defined as the depth direction.
  • Acquisition of the measurement signal by scanning the observation point in the depth direction of the sample Sm is called A-scan.
  • A-scan is realized by using a wavelength sweep light source.
  • the B-scan refers to a scan in a direction perpendicular to the depth direction of the sample Sm.
  • the OCT system 1 includes a light source 102, a coupler 104, a polarization delay unit 110, a circulator 120, a probe 128, a reference arm 130, a polarization separation detection unit 150, a photodetector 190, and a measurement signal processing device 200.
  • the light source 102, the coupler 104, the polarization delay unit 110, the circulator 120, the probe 128, the polarization separation detection unit 150, and the photodetector 190 are components that constitute the optical system, respectively.
  • the components of the optical system are connected using an optical fiber as an optical path.
  • the light source 102 is a wavelength sweep light source (Wavelength Swept Source) that periodically generates light having a wavelength that is swept within a predetermined wavelength width (for example, 40 to 120 nm).
  • the light source 102 has a near-infrared wavelength (for example, 1000 to 1400 nm) such as SLD (Superluminescent Diode).
  • SLD Superluminescent Diode
  • the coupler 104 separates the light incident from the light source 102 into two systems, a probe arm and a reference arm 130, at a predetermined intensity ratio.
  • the ratio of the light intensity to the probe arm to the light intensity to the reference arm 130 is, for example, 90%: 10%.
  • the light supplied to the probe arm is supplied to the polarization separation detection unit (PDDU: Polarization Diversity Detection Unit) 150.
  • PDDU Polarization Diversity Detection Unit
  • a fiber collimator 106, a polarization controller 108, a polarization delay unit (PDU: Polarization Delay Unit) 110, a circulator 120, a fiber collimator 122, a polarization controller 124, an objective lens 126, and a probe 128 are connected in that order. Is the route.
  • the probe arm is also called a sample arm or a measurement arm.
  • the light supplied to the probe arm is incident on the polarization delay unit 110 via the fiber collimator 106 and the polarization controller 108.
  • the light of the other system is incident on the PPDU 150 via the reference arm 130.
  • the polarization controller 108 amplifies the intensity of the incident light to a predetermined sufficient intensity, and emits the amplified light.
  • the PDU 110 includes a linear polarizing device (Linear Polarizer) 112, a polarizing beam splitter (PBS: Polarization Beam Splitter) 114, and two right-angled prisms (RAP: Right Angle Prism) 116 and 118.
  • the PDU 110 separates the incident light into a horizontally polarized light component and a vertically polarized light component as components having two polarized states orthogonal to each other, and supplies the light obtained by combining the separated components to the circulator 120.
  • the linear polarizing device 112 converts the polarized state of the light incident from the coupler 104 into linearly polarized light, and emits the converted light to the PBS 114.
  • the polarization angle of the linear polarizing device 112 is set to 45 °.
  • the PBS 114 has a reflective layer whose surface is arranged so that the incident angle is 45 °, and the vertically polarized light component of the incident light incident on the reflective layer from the linear polarizing device 112 is transmitted as transmitted light and reflected.
  • the horizontally polarized light component reflects off the surface of the layer as reflected light.
  • the reflected light containing the horizontally polarized light component and the transmitted light containing the vertically polarized light component from the PBS 114 are incident on the RAPs 116 and 118, respectively.
  • the reflective layer of PBS 114 combines the transmitted light transmitted from the incident light including the horizontally polarized light component incident on the reflective layer from RAP116 and the reflected light with respect to the incident light including the vertically polarized light component incident from PBS 118.
  • the combined light is emitted to the circulator 120.
  • Each of the RAPs 116 and 118 has a shape in which the cross section parallel to the optical path is a right-angled isosceles triangle, and the bases of the right-angled isosceles triangles are arranged in a direction orthogonal to the optical path from the PBS 114.
  • Light incident from PBS 114 passes through a side parallel to the base of a right-angled isosceles triangle, is reflected by one of the two sides facing the base, and the reflected light is parallel to the other. It is reflected by the side surface and returns to the side surface parallel to the bottom surface, passes through the side surface, and is incident on the PBS 114.
  • the position of the RAP 118 is adjusted in advance so that the optical path lengths of the PBS 114 and the RAP 116 and the optical path lengths of the PBS 114 and the RAP 118 are significantly different. As a result, the horizontally polarized light component and the vertically polarized light component incident on the sample Sm from the PBS 114 are superimposed and emitted with a predetermined phase difference between them.
  • the circulator 120 emits the light incident from the PDU 110 to the objective lens 126 via the fiber collimator 122 and the polarization controller 124.
  • the objective lens 126 collects the light incident on its own portion and irradiates the sample Sm via the probe 128.
  • the light reflected, scattered, or both obtained from the sample Sm is converted into parallel light by the objective lens 126 via the probe 128 and converted into parallel light by the objective lens 126 as a measurement beam to the circulator 120 via the polarization controller 124 and the fiber collimator 122. It returns and is incident on the PPDU150.
  • the reference arm 130 is a path formed by connecting a fiber collimator 132, a fiber Bragg grating (FBG) 134, a fiber collimator 136, a delay line 138, a fiber collimator 140, and a polarization controller 142 in that order. Is.
  • FBG fiber Bragg grating
  • the FBG 134 reflects a component of a specific wavelength of the incident light as reflected light, passes through the remaining components, and is incident on the delay line 138 as transmitted light via the fiber collimator 136.
  • the reflected light from the FBG 134 returns to the coupler 104 via the fiber collimator 132, and is incident on the photodetector 190 from the coupler 104.
  • the band of the component reflected from the FBG 134 is sufficiently narrower than the band of light generated by the light source 102.
  • the photodetector 190 detects the intensity of the reflected light from the FBG 136 and outputs the intensity signal indicating the detected intensity to the measurement signal processing device 200 as a trigger signal.
  • the trigger signal is used as a trigger for the A-scan.
  • the wavelength of the light generated by the light source 102 changes periodically within a predetermined wavelength width range, but the timing at which the wavelength reaches the predetermined wavelength is detected by the photodetector 190, and the A-scan is optical at that timing. It is reset by the system control unit 212.
  • the lower limit of the wavelength width of the photodetector 190 is set in advance as the wavelength to be detected. This is because in SS-OCT, the depth of the observation point to be observed is determined by the wavelength of the probe light.
  • the delay line 138 delays the incident light incident from the FBG 134, and emits the delayed light to the PPDU 150 via the fiber collimator 140 and the polarization controller 142.
  • the delay line 138 makes the delay amount with respect to the incident light variable, and adjusts the delay amount so that the optical path length of the probe arm and the optical path length of the reference arm 130 are equal to each other.
  • the polarization controller 142 adjusts the intensity of the incident light to a predetermined intensity, and emits the light with the adjusted intensity.
  • the PPDU150 is a linear polarization device 152, a non-polarization beam splitter (NPB) 154, two PBS 156, 158, four receivers 162, 164, 166, 168 and two balanced polarization detectors.
  • BPD Balanced Polarization Detector
  • the linear polarizing device 152 converts the polarization state of the light incident from the reference arm 130 into linearly polarized light, and emits the converted light to NPBS154.
  • the polarization angle of the linear polarizing device 152 is set to 45 °.
  • the NPBS 154 combines the incident light incident from the reference arm 130 via the linear polarizing unit 152 with the incident light incident from the probe arm.
  • the NPBS 154 has a reflective layer whose surface is arranged in a direction in which the incident angle is 45 ° with respect to each of the incident light from the reference arm 130 and the incident light from the probe arm.
  • the reflective layer combines the transmitted light obtained by transmitting the incident light from the reference arm 130 and the reflected light obtained by reflecting the incident light from the probe arm, and the interference light obtained by the combined wave is combined with the PBS 158.
  • the reflective layer combines the reflected light obtained by transmitting the incident light from the reference arm 130 and the transmitted light obtained by transmitting the incident light from the probe arm, and the interference light obtained by the combined wave is combined with the PBS 156.
  • the PBS 156 separates the interference light incident from the NPBS 154 into a horizontally polarized wave component and a vertically polarized wave component, and emits the separated horizontally polarized wave component and the vertically polarized wave component to the receivers 162 and 166, respectively.
  • the receivers 162 and 166 receive the horizontally polarized wave component and the vertically polarized wave component incident from the PBS 156, respectively, and guide the light guides to the BPD 170 and 172 as the first horizontally polarized wave component and the first vertically polarized wave component, respectively.
  • the first horizontally polarized component and the first vertically polarized component correspond to the horizontally polarized component and the vertically polarized component of the interference light based on the horizontally polarized component incident on the sample Sm, respectively.
  • PBS158 separates the light incident from NPBS154 into a horizontally polarized wave component and a vertically polarized wave component, and emits the separated horizontally polarized wave component and vertically polarized wave component to the receivers 164 and 168, respectively.
  • the receivers 164 and 168 receive the horizontally polarized wave component and the vertically polarized wave component incident from the PBS 158, respectively, and guide the light guides to the BPD 170 and 172 as the second horizontally polarized wave component and the second vertically polarized wave component, respectively.
  • the second horizontally polarized light component and the second vertically polarized light component correspond to the horizontally polarized light component and the vertically polarized light component of the interference light based on the vertically polarized light component incident on the sample Sm, respectively.
  • the BPD 170 detects the first horizontally polarized light component and the second horizontally polarized light component guided from the receivers 162 and 166, respectively, and indicates the intensities of the detected first horizontally polarized light component and the second horizontally polarized light component. It is converted into a first horizontally polarized light spectrum interference signal and a second horizontally polarized light spectrum interference signal which are analog electric signals.
  • the BPD 170 passes the generated first horizontally polarized spectrum interference signal and second horizontally polarized spectrum interference signal via a low-pass filter (LPF: Low Pass Filter) 182 and a high-pass filter (HPF: High Pass Filter) 186. And output to the measurement signal processing device 200.
  • LPF Low Pass Filter
  • HPF High Pass Filter
  • the BPD 172 detects the first vertically polarized signal and the second vertically polarized component guided from the receivers 164 and 168, respectively, and indicates the intensity of the detected first vertically polarized component and the second vertically polarized component. It is converted into a first vertically polarized light spectrum interference signal and a second vertically polarized light spectrum interference signal which are analog electric signals.
  • the BPD172 outputs the generated first vertical polarization spectrum interference signal and second vertical polarization spectrum interference signal to the measurement signal processing device 200 via LPF184 and HPF188.
  • the first horizontally polarized spectrum interference signal, the second horizontally polarized spectrum interference signal, the first vertically polarized spectrum interference signal, and the second vertically polarized spectrum interference signal are used to generate an OCT image of one frame each.
  • Jones matrix tomography JM-OCT: Jones Matrix-OCT
  • the width of the wavelength of the light emitted by the light source 102 is set by the BPD 170 and 172, and the signal values of a predetermined number of samples (for example, 400 to 2000 samples) are sampled at a predetermined sampling frequency for each A-scan.
  • LPF182, 184 and HPF186, 188 are, for example, Chebyshev filters, respectively.
  • the cutoff frequency of LPF182 and 184 is, for example, 62 MHz.
  • the cutoff frequency of HPF186 and 188 is, for example, 1 MHz.
  • FIG. 10 is a block diagram showing a configuration example of the measurement signal processing device 200 according to the present embodiment.
  • the measurement signal processing device 200 is in the sample from the first horizontally polarized spectrum interference signal, the second horizontally polarized spectrum interference signal, the first vertically polarized spectrum interference signal, and the second vertically polarized spectrum interference signal input from the PPDU 150.
  • the polarization characteristics at the observation points of are analyzed, the time-changing characteristics of the polarization characteristic values indicating the analyzed polarization characteristics are analyzed, and the time-variation characteristic values indicating the analyzed time-variation characteristics are determined.
  • the measurement signal processing device 200 converts the signal value into a signal value indicating a color or gradation corresponding to the time fluctuation characteristic value determined for each observation point, and generates image data for each pixel corresponding to the observation point. , The generated image data may be output.
  • the measurement signal processing device 200 includes a control unit 210, a storage unit 230, an input / output unit 240, a display unit 250, and an operation unit 260.
  • a part or all the functions of the control unit 210 are realized as a computer including a processor such as a CPU (Central Processing Unit), for example.
  • the processor reads a program stored in the storage unit 230 in advance, performs a process instructed by a command described in the read program, and performs its function. In the present application, performing a process instructed by a command described in a program may be referred to as executing a program, executing a program, or the like.
  • control unit 210 is not limited to general-purpose hardware such as a processor, and may be configured to include dedicated hardware such as LSI (Large Scale Integration) and ASIC (Application Specific Integrated Circuit).
  • the control unit 210 includes an optical system control unit 212, a measurement signal acquisition unit 214, an ellipsometry unit 216, a fluctuation characteristic analysis unit 218, an image processing unit 220, and an output processing unit 222.
  • the measurement signal acquisition unit 214 and the ellipsometry unit 216 according to the present embodiment function as a measurement unit that acquires an OCT signal and acquires a signal value based on the acquired OCT signal for each observation point of the sample.
  • the fluctuation characteristic analysis unit 218 functions as an evaluation unit for calculating a time fluctuation characteristic value indicating the time fluctuation characteristic of the acquired signal value for each predetermined period.
  • the optical system control unit 212 drives a drive mechanism that changes the position of the probe and scans the observation point of the sample Sm (B-scan).
  • the optical system control unit 212 sets the observation point of the sample Sm in a predetermined direction (for example, the x direction on the front surface orthogonal to the depth direction of the sample Sm) intersecting the depth direction of the sample Sm (hereinafter, z direction). It is scanned.
  • a predetermined direction for example, the x direction on the front surface orthogonal to the depth direction of the sample Sm
  • z direction intersecting the depth direction of the sample Sm
  • the optical system control unit 212 returns the position of the probe 128 to the reference point every time the number of observation points in the x direction from the reference point reaches a predetermined number (number of lines).
  • the measurement signal acquisition unit 214 acquires the measurement signal of the next frame. By repeating the measurement signal, the measurement signal is accumulated for each frame with the passage of time.
  • the x-direction, the z-direction, and the time t are shown in the right, lower, and upper right directions, respectively, and each frame is shown as an individual rectangle.
  • the measurement signal acquisition unit 214 receives a first horizontal polarization spectrum interference signal, a second horizontal polarization spectrum interference signal, a first horizontal polarization spectrum interference signal, and a second vertical polarization from the PPDU 150 via the input / output unit 240.
  • the spectral interference signal is input as a measurement signal.
  • the measurement signal acquisition unit 214 Fourier-converts the first horizontally polarized spectrum interference signal, the second horizontally polarized spectrum interference signal, the first vertically polarized spectrum interference signal, and the second vertically polarized spectrum interference signal, and each of them
  • the first horizontally polarized OCT signal, the second horizontally polarized OCT signal, the first vertically polarized OCT signal, and the second vertically polarized OCT signal indicating complex amplitude are calculated for each observation point.
  • the measurement signal acquisition unit 214 outputs the calculated first horizontally polarized light OCT signal, second horizontally polarized light OCT signal, first vertically polarized light OCT signal, and second vertically polarized light OCT signal to the ellipsometry unit 216.
  • the ellipsometry unit 216 observes based on the first horizontally polarized OCT signal, the second horizontally polarized OCT signal, the first vertically polarized OCT signal, and the second vertically polarized OCT signal input from the measurement signal acquisition unit 214.
  • a polarization characteristic value indicating the polarization characteristic at a point is calculated for each observation point in the sample Sm.
  • the polarization analysis unit 216 constitutes a Jones matrix showing the polarization characteristics of the light wave related to the measurement, and calculates a predetermined polarization characteristic value as an index value representing the polarization characteristics from the constructed Jones matrix.
  • the Jones matrix is a 2-by-2 matrix showing changes in the polarization state.
  • a first Jones vector and a second Jones vector are arranged in the first column and the second column of the Jones matrix, respectively.
  • the first Jones vector and the second Jones vector are acquired by using incident light having polarization components orthogonal to each other. That is, the first Jones vector is a two-dimensional vector containing complex amplitudes represented by the first horizontally polarized OCT signal and the first vertically polarized OCT signal as elements.
  • the second Jones vector is a two-dimensional vector containing complex amplitudes represented by the second horizontally polarized OCT signal and the second vertically polarized OCT signal as elements.
  • Each column of the Jones matrix corresponds to the polarization state related to the incident, and each row corresponds to the detected polarization state.
  • the polarization analysis unit 216 outputs the calculated polarization characteristic value for each observation point to the fluctuation characteristic analysis unit 218.
  • the Jones matrix directly composed of the OCT signal or the conversion coefficient in the frequency domain thereof at this stage may be referred to as a measured Jones matrix.
  • the variation characteristic analysis unit 218 calculates a predetermined time variation characteristic value as an index value indicating the time change characteristic in the observation period which is a preset predetermined period for the polarization characteristic value input from the polarization analysis unit 216.
  • the observation period is typically about 150-600 frames, for example, when the frame rate is 60 frames / sec.
  • the fluctuation characteristic analysis unit 218 outputs the calculated time fluctuation characteristic value to the image processing unit 220. An example of the polarization characteristic value and the time fluctuation characteristic value will be described later.
  • the image processing unit 220 converts the time fluctuation characteristic value input from the fluctuation characteristic analysis unit 218 into a pixel value within a predetermined value range that can be expressed by a bit depth for each pixel using a predetermined function.
  • the image processing unit 220 calculates a pixel value by, for example, adding a predetermined offset value to a multiplication value obtained by multiplying a function value of a sigmoid function with respect to a time fluctuation characteristic value by a predetermined multiple.
  • the function for converting the time fluctuation characteristic value to the pixel value is not limited to the sigmoid function, and the function value obtained for the input value such as a linear function and a logarithmic function is monotonous with respect to the increase of the input value.
  • the image processing unit 220 generates output image data indicating the pixel values converted for each observation point, and outputs the generated output image data to the display unit 250.
  • the image processing unit 220 may store the output image data in the storage unit 230 according to the control signal input from the output processing unit 222.
  • the output processing unit 222 controls the generation or output of output image data indicating the display image based on the operation signal input from the operation unit 260.
  • the operation signal indicates, for example, the necessity of displaying or storing the displayed image, the observation target area, and the like as parameters.
  • the output processing unit 222 may display a setting screen for guiding parameters that can be set by operation, parameter setting, and parameters that can be set on the display unit, and may configure a user interface related to image display. ..
  • the output processing unit 222 when an operation signal indicating the necessity of displaying the display image is input, the output processing unit 222 outputs a control signal indicating the necessity of the display to the image processing unit 220.
  • the image processing unit 220 outputs output image data to the display unit when a control signal indicating display necessity is input from the output processing unit 222, and outputs an output image when a control signal indicating display / rejection is input from the output processing unit 222. Do not output data to the display.
  • the output processing unit 222 When an operation signal indicating the observation target area is input, the output processing unit 222 outputs a control signal indicating the range of the x-coordinate or the y-coordinate in the observation target area to the measurement signal acquisition unit 214.
  • the measurement signal acquisition unit 214 executes the B-scan within the range indicated by the control signal input from the output processing unit 222.
  • the range of the surface of the sample Sm is defined as a parameter.
  • the storage unit 230 stores various data used for processing executed by the control unit 210 and various data acquired by the control unit 210.
  • the storage unit 230 includes, for example, a non-volatile (non-temporary) storage medium such as a ROM (Read Only Memory), a flash memory, and an HDD (Hard Disk Drive).
  • the storage unit 230 includes, for example, a volatile storage medium such as a RAM (RandomAccessMemory) and a register.
  • the input / output unit 240 wirelessly or wiredly connects various data to other devices so that they can input / output.
  • the input / output unit 240 includes, for example, an input / output interface.
  • the input / output unit 240 is connected to, for example, the polarization separation detection unit 150 and the photodetector 190.
  • the display unit 250 displays an image based on the output image data input from the control unit 210.
  • the display unit 250 may include, for example, a liquid crystal display, an organic electroluminescence display, or the like.
  • the operation unit 260 may be configured to include, for example, buttons, knobs, dials, mice, joysticks, and other members that accept user operations and generate operation signals according to the received operations.
  • the operation unit 260 outputs the acquired operation signal to the control unit 210.
  • the operation unit 260 may be an input interface that receives an operation signal wirelessly or by wire from another device (for example, a portable device such as a remote controller).
  • FIG. 11 is a flowchart showing an example of OCT signal processing according to the present embodiment.
  • the measurement signal acquisition unit 214 receives a first horizontally polarized spectrum interference signal, a second horizontally polarized spectrum interference signal, a first vertically polarized spectrum interference signal, and a second vertically polarized spectrum interference signal from the optical system.
  • the first horizontally polarized OCT signal, the second horizontally polarized OCT signal, the first vertically polarized OCT signal, and the second vertically polarized OCT signal are calculated from each of the acquired signals.
  • the ellipsometry unit 216 has a first horizontally polarized wave OCT signal, a second horizontally polarized wave OCT signal, a first vertically polarized wave OCT signal, and a second vertically polarized wave OCT signal acquired by the measurement signal acquisition unit 214.
  • a Jones matrix is constructed for each observation point based on.
  • the polarization analysis unit 216 calculates a predetermined polarization characteristic value from the constructed Jones matrix.
  • Step S106 The variation characteristic analysis unit 218 calculates a predetermined time variation characteristic value from the change characteristic value calculated by the polarization analysis unit 216 in the preset observation period.
  • Step S108 The image processing unit 220 converts the time fluctuation characteristic value for each observation point calculated by the fluctuation characteristic analysis unit 218 into a pixel value in the pixel corresponding to each observation point.
  • Step S110 The image processing unit 220 outputs output image data indicating the converted pixel value to the display unit 250 (image display). Therefore, the display unit 250 visualizes the distribution of the time-varying characteristic value indicating the polarization characteristic in the observation target region in the observation target region.
  • the polarization analysis unit 216 calculates, for example, a phase retardation as a polarization characteristic value.
  • Phase delay is the phase difference between normal and anomalous rays caused by birefringence. That is, the normal light beam and the abnormal light ray have polarization directions orthogonal to each other with respect to the optical axis of the sample, and pass through the sample at different speeds.
  • the phase delay accumulated from the surface to the observation point of interest in the sample is called cumulative phase delay (CPR).
  • Ellipsometry 216 (3) Measurement Jones matrix at the observation point of the depth z in the inverse matrix J m0 -1 Jones matrix in position on the sample surface corresponding to the observation points of the depth z as shown in equation J multiplied by mz calculated cumulative Jones matrix at the observation point (cumulative Jones matrix) J cz.
  • the cumulative Jones matrix Jcz shows the change in the polarization state from the sample surface to the observation point at the depth z.
  • the polarization analysis unit 216 may calculate a local phase delay (LPR: Local Phase Retardation), which is a polarization phase delay in a local depth region. More specifically, the polarization analyzer 216, (4) cumulative Jones matrix at the observation point depth z 2 to the inverse matrix J CZ1 -1 cumulative Jones matrices at the observation point depth z 1 as shown in equation To calculate the local Jones matrix J l12 between the observation points at the depths z 1 and z 2 obtained by multiplying J cz 2. The local Jones matrix J l12 shows the change in the polarization state from the observation point at depth z 1 to the observation point at depth z 2.
  • LPR Local Phase Retardation
  • the thickness ⁇ z 12 from the observation point at depth z 1 to the observation point at depth z 2 may be, for example, the thickness of the tissue to be observed indicated by the operation signal from the operation unit 260, or in the z direction. It may be the distance between observation points.
  • the ellipsometry unit 216 can calculate the phase difference arg ( ⁇ l1 ⁇ l2 * ) of the two eigenvalues ⁇ l1 and ⁇ l2 of the local Jones matrix J l12 as LPR.
  • the polarization analyzer 216 measures Jones matrix J MZ2 depth z 1 multiplied by the observation point depth z 2 to the inverse matrix J MZ1 -1 measurement Jones matrix at the observation point depth z 1, z 2
  • the local Jones matrix J l12 between the observation points of may be calculated.
  • the ellipsometry unit 216 divides the LPR by 2k 0 ⁇ z 12 and birefringence between the observation points at the depths z 1 and z 2, as shown in equation (5).
  • the rate b 12 may be calculated.
  • k 0 indicates the central wavelength of the incident light.
  • the fluctuation characteristic analysis unit 218 calculates the variance of the polarization characteristic value as an example of the time fluctuation characteristic value within a predetermined observation period.
  • the calculated time variation characteristic value indicates the magnitude of the time variation of the polarization state indicated by the polarization characteristic value.
  • ⁇ l 2 (x, z) indicates the variance of LPR.
  • N, ⁇ (x, z, t i), ⁇ (x, z)> indicates each frame number in a predetermined observation period, LPR, the time average value of the LPR.
  • (X, z) indicates the coordinates of the observation point in the xz plane.
  • t i indicates the i-th sample time.
  • Equation (6) exemplifies the variance of LPR, but the variation characteristic analysis unit 218 may calculate the variance of CPR or the variance of birefringence instead of LPR. The variation characteristic analysis unit 218 may calculate the square root of these variances as the standard deviation.
  • the variation characteristic analysis unit 218 may calculate the logarithmic variance or standard deviation of the polarization characteristic value as another example of the time variation characteristic value.
  • the logarithmic value of the time-averaged value is subtracted from the logarithmic value of the polarization characteristic value at each time, so that the constant that is potentially multiplied by the polarization characteristic value is eliminated. Therefore, it can be applied to the evaluation of substantial time fluctuation characteristics.
  • log ⁇ l 2 (x, z) indicates the logarithmic variance of LPR.
  • Equation (7) takes the variance of the logarithmic value of LPR as an example, but the variation characteristic analysis unit 218 uses the logarithmic variance of CPR or the logarithmic variance of the birefringence index instead of the logarithmic variance of LPR. May be calculated.
  • the variation characteristic analysis unit 218 may calculate the square root of these variances as the standard deviation.
  • the fluctuation characteristic analysis unit 218 may calculate, for example, the dynamic contrast of the polarization characteristic value as the time fluctuation characteristic value within a predetermined observation period.
  • the dynamic contrast of the polarization characteristic value corresponds to a value normalized by dividing the standard deviation with respect to the polarization characteristic value by the time average value.
  • ⁇ d indicates the dynamic contrast of LPR.
  • Equation (8) uses the dynamic contrast of LPR as an example, but the fluctuation characteristic analysis unit 218 may calculate the dynamic contrast of CPR or birefringence instead of LPR.
  • the standard deviation is normalized by the time average value, so that the dynamic contrast can be applied to the evaluation of the time variation characteristics that are more substantial than the standard deviation.
  • the variation characteristic analysis unit 218 may calculate the temporal polarization uniformity (TPU) as another example of the time variation characteristic value.
  • TPU temporal polarization uniformity
  • the ellipsometry unit 216 uses the first Jones vector J 1 and the second Jones vector J 2 forming a subspace of the measurement Jones matrix Jmz related to the observation point as the first Stokes, respectively. It is converted into a vector S 1 and a second Stokes vector S 2 (FIG. 18, step S122).
  • the first Jones vector J 1 and the second Jones vector J 2 are configured to include the elements of the first column and the elements of the second column of the Jones matrix Jmz, respectively.
  • the first Stokes vector S 1 and the second Stokes vector S 2 are four-dimensional vectors representing the polarization states represented by the first Jones vector J 1 and the second Jones vector J 2, respectively. It is composed of the element values of the first Jones vector J 1 and the second Jones vector J 2, respectively.
  • the four element values s 10 to s 13 and s 20 to s 23 forming each of the first Stokes vector S 1 and the second Stokes vector S 2 are set to the 0th to third Stokes parameters s 10 to. They are called s 13 , s 20 to s 23.
  • the 0th Stokes parameters s 10 and s 20 are the power of the horizontal component
  • the first Stokes parameter is the difference between the power of the horizontal component
  • 2 that is, between the components orthogonal to each other. Shows the difference.
  • the second Stokes parameter is twice the real part of the complex conjugate product g 1H g 1V * , g 2H g 2V * of the horizontal component and the vertical component, that is, the horizontal components g 1H , g 2H and the vertical component.
  • the third Stokes parameter is twice the imaginary part of the complex conjugate of the horizontal component and the vertical component g 1H g 1V * , g 2H g 2V * , that is, the horizontal components g 1H , g 2H and the vertical component.
  • the fluctuation characteristic analysis unit 218 calculates the time average values ⁇ S 1 > and ⁇ S 2 > within the observation period for each of the first Stokes vector S 1 and the second Stokes vector S 2 , and uses the respective element values.
  • the sum of the 0th to 3rd Stokes parameters is the 0th time average value ⁇ s 10 + s 20 >, the 1st time average value ⁇ s 11 + s 21 >, the 2nd time average value ⁇ s 12 + s 22 >, and the 3rd time. It is calculated as an average value ⁇ s 13 + s 23 > (FIG. 18, step S124).
  • the fluctuation characteristic analysis unit 218 divides the square root of the sum of squares with respect to the first time average value, the second time average value, and the third time average value by the 0th time average value as shown in the equation (11).
  • the obtained value is defined as TPU (FIG. 18, step S126).
  • the TPU takes a smaller value as the time variation of the polarized state is smaller. Therefore, it is expected that the TPU becomes smaller as the activity of the biological tissue used as a sample becomes more active.
  • the above-mentioned time-varying characteristic value for CPR, LPR or birefringence becomes larger as the time-varying of the polarized state becomes remarkable.
  • the measurement Jones matrix contains a noise component as shown in Eq. (12).
  • n 1H , n 1V , n 2H , and n 2V indicate noise components added to the signal components E 1H , E 1V , E 2H , and E 2V, respectively.
  • the ellipsometry unit 216 has the time average power of the element values included in the time average value ⁇ S 1> of the first Stokes vector and the time average value ⁇ S> of the second Stokes vector.
  • the noise component can be compensated by subtracting the time average power of the element value of the noise component from.
  • the fluctuation characteristic analysis unit 218 replaces the time average value ⁇ S 1 > of the first Stokes vector and the time average value ⁇ S 2 > of the second Stokes vector before correction of the noise component in the equation (11).
  • the TPU from which the noise component has been removed is calculated by substituting the time average value ⁇ S 1 '> of the first Stokes vector and the time average value ⁇ S 2 '> of the second Stokes vector with respect to the signal component after correction of the noise component. be able to.
  • the fluctuation characteristic analysis unit 218 may calculate, for example, the von Neumann Entropy of the Jones matrix determined at the observation point as the time fluctuation characteristic value within the predetermined observation period.
  • the von Neumann entropy is defined by Eq. (14).
  • ⁇ i indicates the eigenvalues of the Hermitian matrix T (described later), and the individual eigenvalues ⁇ i are normalized by the sum of those eigenvalues to obtain the normalized eigenvalues ⁇ i '.
  • LPR expected value E (R) as shown in equation is known to be a weighted average value of R i is a LPR in the optical axis corresponding to an individual eigenvalues. Weighting factors to be multiplied to the LPR of the individual optical axes of the weighted average is given by the normalized eigenvalue lambda i 'corresponding to the optical axis. Therefore, it can be seen as a kind of time varying characteristic value Von Neumann entropy of the Jones matrix defined with normalized eigenvalue lambda i 'also according to the polarization state. The von Neumann entropy of the Jones matrix is described in detail in the following literature, but in the present embodiment, the randomness due to the time variation of the Jones matrix is evaluated as the polarization characteristic value.
  • the polarization analysis unit 216 and the variation characteristic analysis unit 218 can calculate the von Neumann entropy H by the following procedure.
  • the ellipsometry unit 216 is a four-dimensional vector [g 1H g 2H g 1V g 2V ] formed by rearranging the elements of the measurement Jones matrix Jmz in the order of each row and column for each observation point.
  • T is configured as the target vector ⁇ L.
  • This target vector kappa L shows substantially the same value as the measurement Jones matrix Jmz.
  • the ellipsometry unit 216 multiplies the Hermitian conjugate ⁇ L + of the target vector by the target vector ⁇ L to calculate a 4-by-4 square matrix ⁇ L ⁇ L +.
  • the variation characteristic analysis unit 218 defines the time average value ⁇ L ⁇ L + > of the square matrix ⁇ L ⁇ L + within a predetermined period as the matrix T.
  • the matrix T is a 4-by-4 semi-fixed Hermitian matrix (Positive Semidefinite Hermitian Matrix).
  • the von Neumann entropy has a value of 0 or more and 1 or less.
  • the von Neumann entropy H is 1.
  • the base of the logarithmic value is 4.
  • the fluctuation characteristic analysis unit 218 subtracts the von Neumann entropy H n of the noise component from the von Neumann entropy H m determined based on the measurement Jones matrix, and subtracts the von Neumann entropy H n of the signal component from the von Neumann entropy H m of the signal component.
  • Entropy H s may be defined.
  • Variation characteristic analyzer 218, (18) as shown in equation Von Neumann entropy H n (E 1, E 2) first Jones vector E 1 which forms a subspace of the Jones matrix J mz a noise component, a second Of the first noise component entropy H (E 1 ), which is the entropy of the first noise component n 1 applied to each of the Jones vectors E 2 , and the second noise component entropy H (E 2 ), which is the entropy of the second noise component n 2. It can be approximated to the sum. However, it is assumed that the first noise component n 1 and the second noise component n 2 are independent of each other.
  • the fluctuation characteristic analysis unit 218 sets the i-th noise component entropy H (E i ) as the j-eigenvalue ⁇ j (i) of the i-noise component and its logarithmic log ( ⁇ j ( ⁇ j (). i) The sum of) can be calculated by reversing the positive and negative signs.
  • the fluctuation characteristic analysis unit 218 adds or subtracts P (i), which is the TPU of the i-Jones vector, to 1, and divides the value by 2 to the first eigenvalue ⁇ . It can be calculated as 1 (i) and the second eigenvalue ⁇ 2 (i).
  • the fluctuation characteristic analysis unit 218 corrects P (i) , which is the TPU of the i-Jones vector, with the above noise component, and then corrects the first Stokes vector and the second Stokes vector, respectively. It can be calculated by dividing the square root of the sum of squares of the time average values of the first to third Stokes parameters by the time average value of the 0th Stokes parameter.
  • FIG. 13 is a diagram showing an example of the spatial distribution of the dispersion of the birefringence index.
  • the dispersion of the birefringence index for each observation point in the living tissue is shown by shading.
  • the observation points are distributed in the xz plane.
  • the upper part filled with black indicates the outside of the tissue. As a whole, the dispersion of the birefringence tends to be larger inside than on the surface of the structure.
  • FIG. 14 is a diagram showing another example of the spatial distribution of the dispersion of the birefringence index.
  • FIG. 14 shows the dispersion of the birefringence index for each observation point in the living tissue in shades.
  • FIG. 15 shows the mean local birefringence for each observation point.
  • the same biological tissue as in FIG. 14 is the observation target.
  • a biological tissue different from that in FIG. 13 is the observation target.
  • the portion brighter than the surroundings shows the distribution range of the living tissue.
  • FIG. 16 is a diagram showing an example of the correlation between the dispersion of the birefringence index and the average local birefringence index.
  • the vertical axis and the horizontal axis of FIG. 16 show the variance of the birefringence and the average local birefringence, respectively.
  • FIG. 16 shows that there is a significant correlation between the local birefringence index and the variance of the birefringence index. The correlation coefficient was 0.776.
  • FIG. 17 is a diagram showing an example of the correlation between the dispersion of the birefringence index and the dispersion of the logarithmic intensity.
  • the vertical axis and the horizontal axis of FIG. 17 show the dispersion of the birefringence index and the dispersion of the logarithmic intensity, respectively.
  • the logarithmic strength is the logarithmic value of the signal strength for each observation point.
  • FIG. 17 shows that there is no significant correlation between the local birefringence index and the variance of the birefringence index.
  • the correlation coefficient was 0.280. This confirms that the polarized state cannot be fully expressed by the logarithmic intensity alone.
  • FIG. 19 is a diagram showing an example of the spatial distribution of TPU.
  • FIG. 19 shows the TPU for each observation point in the living tissue in shades.
  • the same biological tissue as in FIG. 13 is the observation target.
  • the TPU tends to be smaller inside than on the surface of the tissue. This tendency is opposite to the dispersion of the birefringence index. This is because the larger the TPU, the smaller the time change of the polarized state, whereas the larger the dispersion of the birefringence, the larger the time change of the polarized state.
  • FIG. 20 is a diagram showing an example of the time change of the dispersion of the birefringence index.
  • FIG. 21 is a diagram showing another example of the time change of TPU.
  • FIG. 22 is a diagram showing an example of a time change of the dynamic contrast of the birefringence index.
  • FIG. 23 is a diagram showing an example of time change of logarithmic intensity dispersion. 20, FIG. 21, FIG. 22, and FIG. 23 show hourly birefringence variance, TPU, birefringence index, and logarithmic intensity variance at a certain observation point in a common biological tissue, respectively. The activity state of the living tissue decreases with the passage of time.
  • FIG. 20 shows that the dispersion of the birefringence tends to decrease significantly with the passage of time from time 0 to 44 hours, whereas the birefringence is almost the same from 45 hours to 60 hours. It becomes constant.
  • FIG. 20 shows a tendency that the dispersion of the birefringence index significantly decreases with the passage of time from time 0 to 44 hours later.
  • the correlation coefficient was -0.9486.
  • the birefringence was almost constant from the lapse of 45 hours to the lapse of 60 hours.
  • the correlation coefficient was -0.1711.
  • FIG. 21 shows a tendency that the TPU increases significantly with the passage of time from time 0 to the passage of 42 hours.
  • the correlation coefficient was 0.9413.
  • the TPU was almost constant from the lapse of 43 hours to the lapse of 60 hours.
  • the correlation coefficient was 0.0735.
  • FIG. 22 shows a tendency for the dynamic contrast of the birefringence index to decrease significantly with the passage of time. The correlation coefficient was ⁇ 0.905.
  • the logarithmic intensity variance decreases from time 0 to 3 hours, the logarithmic intensity variance increases from 4 hours to 19 hours, and the logarithm increases from 20 hours to 42 hours.
  • the intensity dispersion decreased, and the logarithmic intensity dispersion became almost constant from the lapse of 42 hours to the lapse of 60 hours.
  • the time-varying characteristic value of the polarization characteristic value can be measured for each observation point in the sample, and the activity state of the living tissue is evaluated based on the measured time-varying characteristic value. Show that you can.
  • the measurement signal processing apparatus 200 reflects or scatters from the sample with respect to the first incident component (for example, the horizontally polarized component) incident on the sample in the first polarized state.
  • the first measurement signal for example, the first horizontally polarized spectrum interference signal
  • the first polarization state for example, horizontally polarized light
  • a second measurement signal having a two-polarized state for example, vertically polarized light
  • a second incident component for example, a vertically polarized component incident on a sample in the second polarized state.
  • the third measurement signal having the first polarization state for example, the first vertically polarized spectrum interference signal
  • the second interference for example, the third measurement signal having the first polarization state
  • a polarization analysis unit 216 that determines a polarization characteristic value based on a polarization characteristic at an observation point in a sample based on a fourth measurement signal having a second polarization state for a component (for example, a second vertical polarization spectrum interference signal), and polarization.
  • a variation characteristic analysis unit 218 for determining a time variation characteristic value indicating a time variation characteristic of the characteristic value is provided. With this configuration, the time-varying characteristic value of the polarization characteristic value is determined for each observation point in the sample. Since the time variation of the polarization characteristic value has a significant correlation with the tissue activity, the tissue activity can be quantitatively evaluated by the distribution of the time variation characteristic value determined for each observation point.
  • the polarization analysis unit 216 determines a Jones matrix for each observation point based on the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal, and the Jones matrix and the sample at the observation points in the sample.
  • the cumulative Jones matrix at the observation point may be determined from the Jones matrix on the surface of the surface, and the CPR, which is the phase difference between the eigenvalues of the cumulative Jones matrix, may be determined as the polarization characteristic value.
  • the polarization analysis unit 216 determines a Jones matrix for each observation point based on the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal, and Jones at the first observation point in the sample. Determine the local Jones matrix between the matrix and the Jones matrix at the second station in the sample to the first and second stations, and determine the LPR, which is the phase difference between the two eigenvalues of the local Jones matrix. May be good.
  • the time variation characteristic value of LPR corresponding to the phase difference between the polarized light components generated in the section between the first observation point and the second observation point can be used for the evaluation of the tissue activity. Therefore, the activity of the tissue can be evaluated for each minute region.
  • the ellipsometry unit 216 may determine the birefringence index by dividing the local phase delay by the wave number of the incident light incident on the sample and the thickness of the first observation point and the second observation point. Therefore, the activity of the tissue can be evaluated for each minute region, and the birefringence can be easily compared with other tissues to be observed or an existing substance by using the birefringence.
  • the variation characteristic analysis unit 218 may calculate the time variation characteristic value based on the variance or standard deviation of the polarization characteristic value. Therefore, the magnitude of the time variation of the polarization characteristic value is quantitatively evaluated.
  • the variation characteristic analysis unit 218 may calculate the time variation characteristic value based on the logarithmic variance of the polarization characteristic value or the standard deviation. In the process of calculating the variance or standard deviation, the constant multiplied by the polarization characteristic value is eliminated, so that the substantial time variation of the polarization characteristic value is evaluated. In addition, by taking logarithmic values, it becomes easy to compare time-varying characteristic values with different scales with other tissues to be observed or existing substances.
  • the fluctuation characteristic analysis unit 218 may calculate the dynamic contrast by dividing the standard deviation of the polarization characteristic value by the average value of the birefringence index. By dividing by the average value of the polarization characteristic values, the standard deviation of the polarization characteristic values is normalized, so that the substantial time variation of the polarization characteristic values is evaluated without changing the scale.
  • the polarization analysis unit 216 uses the first measurement signal and the first Jones vector based on the second measurement signal, and the third measurement signal and the second Jones vector based on the fourth measurement signal as polarization characteristic values, respectively. Converted to the first Stokes vector and the second Stokes vector, the fluctuation characteristic analysis unit 218 calculates the TPU as the time fluctuation characteristic value based on the time average value of the first Stokes vector and the time average value of the second Stokes vector. You may decide. According to this configuration, the uniformity of the polarization state at the observation point over time is quantified by using the TPU. TPU tends to increase as the activity of the tissue decreases. Therefore, the inactive state of the tissue can be quantitatively evaluated by the distribution of TPU determined for each observation point.
  • the fluctuation characteristic analysis unit 218 describes the time average value of the corrected first Stokes vector obtained by subtracting the noise component from the first Stokes vector and the time average of the corrected second Stokes vector obtained by subtracting the noise component from the second Stokes vector.
  • the TPU may be determined based on the value. According to this configuration, the noise component is removed from the first Stokes vector and the second Stokes vector, and the signal component is left. Therefore, the influence of noise on the TPU can be suppressed, and the activity of the tissue can be accurately evaluated.
  • the polarization analysis unit 216 determines a Jones matrix for each observation point from the first measurement signal, the second measurement signal, the third measurement signal, and the fourth measurement signal as the polarization characteristic value, and the variation characteristic analysis unit. 218 may calculate the von Neumann entropy of the Jones matrix as the time fluctuation characteristic value. According to this configuration, the randomness of the Jones matrix indicating the polarization state at the observation point is quantified. Therefore, the tissue activity can be quantitatively evaluated by the distribution of von Neumann entropy determined for each observation point.
  • the fluctuation characteristic analysis unit 218 has converted the first measurement signal and the first Jones vector based on the second measurement signal, and the third measurement signal and the second Jones vector based on the fourth measurement signal, respectively.
  • the entropy of the noise component may be calculated from the time polarization uniformity of the 1st Stokes vector and the time polarization uniformity of the 2nd Stokes vector, and the von Neumann entropy may be corrected based on the entropy of the noise component. According to this configuration, the contribution of the entropy of the noise component is compensated from the von Neumann entropy, and the von Neumann entropy of the signal component is obtained. Therefore, the influence of noise on the von Neumann entropy can be suppressed, and the activity of the tissue can be accurately evaluated.
  • the measurement signal processing device 200 is an image processing unit that generates image data having an output value as a signal value with respect to the time fluctuation characteristic value for each observation point by using a function that gives an output value that changes monotonically with respect to a change in the input value. 220 may be provided. According to this configuration, an image having a brightness or color tone distribution corresponding to the time variation characteristic value for each observation point can be obtained.
  • the measurement signal processing device 200 may include an output processing unit 222 that determines an evaluation value indicating an active state of a biological tissue as a sample based on a time fluctuation characteristic value.
  • the evaluation value may be, for example, a real value that increases as the activity increases.
  • the output processing unit 222 for example, a function indicating the relationship between the evaluation value and the time fluctuation characteristic value and its parameters are set in advance.
  • the output processing unit 222 may store the evaluation value in the storage unit 230, or may output the evaluation value to another device.
  • the image processing unit 220 may convert the evaluation value calculated for each observation point by the output processing unit 222 into a signal value for each pixel as described above, and generate image data having the converted signal value. Therefore, the user can easily evaluate the activity of the tissue in the observation region by coming into contact with the obtained image.
  • the evaluation device 20 and the measurement signal processing device 200 are each a part of the OCT system 1 is taken as an example, but the present invention is not limited to this.
  • the evaluation device 20 and the measurement signal processing device 200 may be a single device that is independent of the OCT system 1 and does not have an optical system.
  • the optical system control unit 212 may be omitted in the control unit 210 of the measurement signal processing device 200.
  • the evaluation device 20 may omit the control means for controlling the optical system.
  • the evaluation device 20 and the measurement signal acquisition unit 214 are not limited to the optical system, but acquire detection signals and measurement signals from other devices such as a data storage device and a PC wirelessly or by wire, for example, via a network. May be good.
  • the measurement signal processing device 200 may include any of the input / output unit 240, the display unit 250, and the operation unit 260 as described above, or a part or all of them may be omitted. Further, the evaluation device 20 may also include any of the functional configurations corresponding to the input / output unit 240, the display unit 250, and the operation unit 260, or a part or all of them may be omitted. Further, in the control unit 210 of the measurement signal processing device 200, one or both of the image processing unit 220 and the output processing unit 222 may be omitted.
  • the control unit 210 wirelessly or wiredly transfers the generated data indicating the time fluctuation characteristic value to another device such as a data storage device, a PC, or another image processing device, for example. It may be output via the network.
  • the evaluation device 20 may also include any of functional configurations corresponding to one or both of the image processing unit 220 and the output processing unit 222, or a part or all of them may be omitted.
  • the device as the output destination has the same function as the image processing unit 220, that is, the output image data is generated based on the data input from the evaluation device 20 or the measurement signal processing device 200, and the image based on the generated output image data. May have a function of displaying.
  • the functional configuration corresponding to one or both of the image processing unit 220 and the output processing unit 222 of the measurement signal processing device 200 may be omitted.
  • a part or all of the evaluation device 20 or the measurement signal processing device 200 in the above-described embodiment may be realized as an integrated circuit such as an LSI (Large Scale Integration).
  • LSI Large Scale Integration
  • Each functional block of the evaluation device 20 or the measurement signal processing device 200 may be individually made into a processor, or a part or all of them may be integrated into a processor.
  • the method of making an integrated circuit is not limited to LSI, and may be realized by a dedicated circuit or a general-purpose processor.

Landscapes

  • Health & Medical Sciences (AREA)
  • Life Sciences & Earth Sciences (AREA)
  • Engineering & Computer Science (AREA)
  • Physics & Mathematics (AREA)
  • Chemical & Material Sciences (AREA)
  • Biomedical Technology (AREA)
  • Pathology (AREA)
  • General Health & Medical Sciences (AREA)
  • General Physics & Mathematics (AREA)
  • Immunology (AREA)
  • Biochemistry (AREA)
  • Analytical Chemistry (AREA)
  • Urology & Nephrology (AREA)
  • Hematology (AREA)
  • Molecular Biology (AREA)
  • Food Science & Technology (AREA)
  • Medicinal Chemistry (AREA)
  • Biophysics (AREA)
  • Optics & Photonics (AREA)
  • Investigating Or Analysing Materials By Optical Means (AREA)

Abstract

生体組織のダイナミクスを定量的に評価することを可能とする生体組織の評価方法を提供する。本実施形態の評価方法は、試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得し、前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する。本実施形態は、評価装置またはプログラムとしても実現することができる。

Description

評価装置、評価方法およびプログラム
 本発明は、評価装置、方法およびプログラムに関する。例えば、光干渉断層計(OCT:Optical Coherence Tomography)で計測された計測信号を処理することにより、生体組織などの試料の状態を可視化ならびに定量評価するための技術に関する。
 本願は、2019年11月15日に日本に出願された特願2019-207348号と、2020年4月9日に日本に出願された特願2020-070309号と、について優先権を主張し、それらの内容をここに援用する。
 近年、OCTを用いてex vivo試料をイメージングする「OCT顕微鏡」と呼ばれる技術が研究されている。しかし、一般にOCT顕微鏡は形態イメージングする技術であり、代謝等の組織の機能をイメージングすることはできない。これに対し、「ダイナミックOCT」と呼ばれる信号解析手法が提案されている(非特許文献1)。しかし、この手法は定量性に乏しく、この手法で生体の活動の度合いを正しく評価することは困難である。また、この手法は、Full-field OCT(FF-OCT)と呼ばれる特殊なタイプのOCTに適合しており、広く使われている走査型OCTでの実装は困難である。
 本発明は上記事情に鑑みてなされたものであり、試料の動的特性、例えば、生体組織のダイナミクスや細胞内活動を定量的に評価することを可能とする評価方法を提供することを課題の一つとする。
 上記課題を解決するため、本発明は以下の手段を採用している。
 (1)本発明の一態様は、試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測部と、前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価部と、を有する評価装置である。
 (2)本発明の他の態様において、前記評価部は、前記時間変動特性値として、前記信号値の分散を算出してもよい。
 (3)本発明の他の態様において、前記評価部は、前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差の二乗和を前記所定期間におけるフレーム数で除算して前記分散を観測点ごとに算出してもよい。
 (4)本発明の他の態様において、前記評価部は、前記信号値と、当該信号値を時間シフト量τで時間シフトした信号値との相関係数を前記時間シフト量τごとに算出し、前記時間シフト量の増加に伴う前記相関係数の減衰速度を前記時間変動特性値として算出してもよい。
 (5)本発明の他の態様において、前記評価部は、前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差の二乗和を分散として算出し、前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差と、当該フレーム時刻から時間シフト量τで時間シフトしたシフト時刻における前記OCT信号の信号強度と、時間シフトした当該信号強度の平均値との偏差との積の総和を共分散として算出し、前記共分散を前記分散で除算して前記相関係数として、前記時間シフト量τごとに算出し、前記時間シフト量τごとの前記相関係数を用いて所定の減衰関数を用いて回帰分析し、前記相関係数を近似する当該減衰関数のパラメータを前記減衰速度として観測点ごとに算出してもよい。
 (6)本発明の他の態様において、前記評価部は、前記時間シフト量τを非零として算出される前記相関係数を用いて、前記減衰速度を算出してもよい。
 (7)本発明の他の態様において、前記計測部は、前記試料に第1偏光状態で入射される第1入射成分に対して前記試料から反射または散乱した成分に前記第1入射成分が干渉した第1干渉成分のうち第1偏光状態を有する第1計測信号、前記第1干渉成分に対する第2偏光状態を有する第2計測信号、前記試料に第2偏光状態で入射される第2入射成分に対して前記試料から反射または散乱した成分に前記第2入射成分が干渉した第2干渉成分のうち第1偏光状態を有する第3計測信号、および前記第2干渉成分に対する第2偏光状態を有する第4計測信号に基づいて、前記試料内の観測点における偏光特性に基づく偏光特性値を定め、前記評価部は、前記偏光特性値の時間変動特性を示す前記時間変動特性値を定めてもよい。
 (8)本発明の他の態様において、前記計測部は、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、前記試料内の観測点におけるジョーンズ行列と前記試料の表面におけるジョーンズ行列から前記観測点における累積ジョーンズ行列を定め、前記偏光特性値として、前記累積ジョーンズ行列の固有値間の位相差である累積位相遅指標値を定めてもよい。
 (9)本発明の他の態様において、前記計測部は、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、前記試料内の第1観測点におけるジョーンズ行列と前記試料内の第2観測点におけるジョーンズ行列から前記第1観測点ならびに前記第2観測点との間の局所ジョーンズ行列を定め、前記局所ジョーンズ行列の固有値間の位相差である局所位相遅延に基づく前記偏光特性値を定めてもよい。
 (10)本発明の他の態様において、前記計測部は、前記局所位相遅延を前記試料に入射する入射光の波数と、前記第1観測点と前記第2観測点との厚みで除算して複屈折率を定めてもよい。
 (11)本発明の他の態様において、前記評価部は、前記偏光特性値の分散または標準偏差に基づく前記時間変動特性値を算出してもよい。
 (12)本発明の他の態様において、前記評価部は、前記偏光特性値の対数値の分散または標準偏差に基づく前記時間変動特性値を算出してもよい。
 (13)本発明の他の態様において、前記評価部は、前記偏光特性値の標準偏差を前記複屈折率の平均値で除算して動的コントラストを算出してもよい。
 (14)本発明の他の態様において、前記計測部は、前記偏光特性値として、前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルを、それぞれ第1ストークスベクトル、および第2ストークスベクトルに変換し、前記評価部は、前記時間変動特性値として、前記第1ストークスベクトルの時間平均値と前記第2ストークスベクトルの時間平均値に基づいて時間偏光均一度を定めてもよい。
 (15)本発明の他の態様において、前記計測部は、前記計測部は、前記第1ストークスベクトルからノイズ成分を差し引いた補正後の第1ストークスベクトルの時間平均値と、前記第2ストークスベクトルからノイズ成分を差し引いた補正後の第2ストークスベクトルの時間平均値に基づいて時間偏光均一度を定めてもよい。
 (16)本発明の他の態様において、前記計測部は、前記偏光特性値として、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、前記評価部は、前記時間変動特性値として、前記ジョーンズ行列のフォンノイマンエントロピーを算出してもよい。
 (17)本発明の他の態様において、前記評価部は、前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルから、それぞれ変換された第1ストークスベクトルの時間偏光均一度、および第2ストークスベクトルの時間偏光均一度からノイズ成分のエントロピーを算出し、前記フォンノイマンエントロピーを前記ノイズ成分のエントロピーに基づいて補正してもよい。
 (18)本発明の他の態様において、前記第1偏光状態は水平偏光であり、前記第2偏光状態は垂直偏光であり、前記第1計測信号は第1水平偏波スペクトル干渉信号であり、前記第2計測信号は第2水平偏波スペクトル干渉信号であり、前記第3計測信号は第1垂直偏波スペクトル干渉信号であり、前記第4計測信号は第2垂直偏波スペクトル干渉信号であってもよい。
 (19)本発明の他の態様において、前記評価部は、前記所定期間よりも長い観測期間間隔ごとに前記時間変動特性値を算出してもよい。
 (20)本発明の他の態様は、前記時間変動特性値に基づいて前記試料の活性状態を示す評価値を定める出力処理部を備えてもよい。
 (21)本発明の他の態様は、入力値の変化に対して単調に変化する出力値を与える関数を用いて前記観測点ごとの前記時間変動特性値に対する出力値を信号値として有する画像データを生成する画像処理部を備えてもよい。
 (22)本発明の他の態様は、評価装置における方法であって、試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測ステップと、前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価ステップと、を有する評価方法である。
 (23)本発明の他の態様は、評価装置のコンピュータに、試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測手順と、前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価手順と、を実行させるためのプログラムである。
 本実施形態によれば、従来は見過ごされていた微小な変化(ゆらぎ)を検出することができ、生体組織のダイナミクスの定量的な評価を実現することができる。
第1の実施形態に係る生体組織の評価方法について説明する図である。 生体組織の一部に対し、比較例1、実施例1、2の方法を適用して得られた画像である。 実施例1の画像に対応するスペックルバリアンス(スペックル信号強度の時間分散)を示す図である。 実施例2の画像に対応するOCT相関の減衰速度を示す図である。 実施例3の方法を適用した場合のOCT信号強度の減衰係数を示す図である。 実施例1の方法を適用した生体組織の生死状態を示す図である。 実施例2の方法を適用した生体組織の生死状態を示す図である。 生体組織の一部に対し、比較例2、実施例4、5の方法を適用して得られた画像である。 実施例4の画像に対応するスペックルバリアンスを示す図である。 実施例5の画像に対応するOCT相関の減衰速度を示す図である。 実施例6の方法を適用した場合のOCT信号強度の減衰係数を示す図である。 実施例4の方法を適用した生体組織の生死状態を示す図である。 実施例5の方法を適用した生体組織の生死状態を示す図である。 第1の実施形態に係るOCTシステムの一例を示す構成図である。 第2の実施形態に係るOCTシステムの一例を示す構成図である。 第2の実施形態に係る計測信号処理装置の構成例を示すブロック図である。 第2の実施形態に係るOCT信号処理の一例を示すフローチャートである。 第2の実施形態に係る観測期間の説明図である。 複屈折率の分散の空間分布の一例を示す図である。 複屈折率の分散の空間分布の他の例を示す図である。 平均局所複屈折率の分散の空間分布の一例を示す図である。 複屈折率の分散と平均局所複屈折率との相関性の一例を示す図である。 複屈折率の分散と対数強度の分散との相関性の一例を示す図である。 第2の実施形態に係るTPUの算出処理の一例を示すフローチャートである。 TPUの空間分布の一例を示す図である。 複屈折率の分散の時間変化の一例を示す図である。 TPUの時間変化の他の例を示す図である。 複屈折率の動的コントラストの時間変化の一例を示す図である。 対数強度分散の時間変化の一例を示す図である。
 以下、本発明を適用した実施形態に係る評価方法について、図面を用いて詳細に説明する。なお、以下の説明で用いる図面は、特徴をわかりやすくするために、便宜上特徴となる部分を拡大して示している場合があり、各構成要素の寸法比率などが実際と同じであるとは限らない。また、以下の説明において例示される材料、寸法等は一例であって、本発明はそれらに限定されるものではなく、その要旨を変更しない範囲で適宜変更して実施することが可能である。
<第1の実施形態>
 図1は、本発明の第1の実施形態に係る生体組織の評価方法について説明する図である。本実施形態の生体組織の評価方法は、OCTシステム1において、所定の時間内に、試料Smに対して複数回のOCT(光干渉断層)撮影を行う撮影手段(撮影部)10と、各回のOCTの画像から得られるOCT信号強度の時間変化を計測する計測手段(計測部)22と、時間変化に基づいて、試料Smとする生体組織の活動を定量評価する評価手段(評価部)24と、を備えた生体組織の評価装置20を用い、次の手順に沿って、生体組織が有している微小な揺らぎを可視化するとともに、定量的に評価するものである(図8参照)。
 まず、撮影手段10は、生体組織の同じ位置において、同じ部位に対し、所定期間内(例えば、1~15分、典型的には3分以内)に複数回(数回から数百回=数フレームから数百フレーム)、好ましくは10フレーム以上、より好ましくは15フレーム以上、OCT(光干渉断層)撮影(計測)を行う。1回のOCT撮影に係る所定期間(観測期間)は、隣接する観測期間までの観測期間間隔よりも短く、かつ、その時点におけるOCT撮影により得られるOCT信号に表れる試料の動的特性の評価に要求される精度を確保できる期間であればよい。より具体的には、試料の動的特性の周波数成分が、観測期間に対応する最低周波数からフレーム間隔に対応する最高周波数までの周波数帯域に含まれるように、フレーム間隔との間で観測期間を定めておけばよい。また、観測期間間隔は、試料とする生体組織の活動状態の大局的な変化傾向の評価に要求される精度を確保できる期間であればよい。例えば、母体から試料とする生体組織が母体から切り離され細胞死に至るまでの一連の過程(アポトーシス(apoptosis))や、母体において試料とする生体組織に原因が発生してから細胞死に要する一連の過程(壊死(necrosis))などの所要時間よりも十分に短い時間を観測期間間隔として定めておけばよい。試料全体の特性を知る観点から、同じOCT撮影を、他の部位に対しても行うことが好ましい。
 次に、計測手段22と、評価手段24は、各回のOCTの画像から得られるOCT信号強度の時間変化に基づいて、生体組織の活動を定量評価する。定量評価の具体的な方法としては、例えば、スペックルバリアンス(SV: Speckle Variance)、OCT相関の減衰速度(OCDS:OCT Correlation Decay Speed)等の算出、ひいては算出結果の評価が挙げられる。
 スペックルバリアンスは、短時間でのOCT信号強度の分散(揺らぎ)を示すものであり、下記(1)式を用いて算出することができる。
Figure JPOXMLDOC01-appb-M000001
 上記(1)式において、x,zは、それぞれ生体組織の表面内、表面から深さ方向の位置を示し、I(x,z,t)は、線形スケールもしくは対数スケール等で表示された各位置、各時刻におけるOCT信号強度を示し、<I>は平均のOCT信号強度を示し、Nは、所定期間内におけるOCT信号のフレーム数を示している。x,zは、OCT画像をなす個々の画素に対応する観測点の位置に相当する。
 より詳細には、まず、計測手段22は、各フレームの時刻t,t,・・・tにおいて、観測点(x,z)ごとに、それぞれOCT信号強度I(x,z,t)、I(x,z,t)、・・・I(x,z,t)を測定する。続いて、評価手段24は、N回分のOCT信号強度の平均<I>を算出する。続いて、評価手段24は、各回のOCT信号強度と平均<I>との差を2乗する。すなわち、評価手段24は、[I(x,z,t)―<I>]、[I(x,z,t)―<I>]、・・・[I(x,z,t)―<I>]を算出する。最後に、評価手段24は、これらの和をNで割ることにより、スペックルバリアンスσ(x,z)を得ることができる。
 OCT相関の減衰速度は、近接する時刻間(t、t+τ)の相関係数が時間シフト量τの増加に応じて小さくなる速度を示すものであり、下記(2)式を用いて相関係数ρ(x,z,τ)を算出することができる。
Figure JPOXMLDOC01-appb-M000002
 上記(2)式において、x,zは、それぞれ生体組織の表面内、表面から深さ方向の位置を示し、σcov (x,z,τ)、σ (x,z)は、それぞれ、OCT信号強度I(x、z、t)とI(x、z、t+τ)との共分散、スペックルバリアンス(分散)を示している。OCT信号強度I(x、z、t+τ)は、OCT信号強度I(x、z、t)を、時間シフト量τをもって時間シフトして得られる信号値である。即ち、相関係数ρ(x,z,τ)は、観測点(x,z)におけるOCT信号の信号値I(x,z,t)の自己相関関数(auto-correlation function)に相当する。
 より詳細には、評価手段24は、近接する時刻間(t,t+τ)での共分散σcov (x,z,τ)として、[I(x,z,t)―<I>]×[I(x,z,t+τ)―<I>]を算出する。続いて、評価手段24は、この共分散σcov (x、z、τ)を上記(1)式で算出したスペックルバリアンスσ(x,z)で割ることにより、相関係数ρ(x,z,τ)を得ることができる。そして、評価手段24は、得られた相関係数ρ(x,z,τ)に対して回帰分析(regression)を行うことにより時間シフト量τの所定の関数に当てはめて、その関数の関数値が相関係数(x,z,τ)により近似するように、その関数のパラメータをOCT相関の減衰速度として得ることができる。所定の関数は、時間シフト量の増加に応じて、その時間シフト量に対して与えられる関数値が減衰する関数、例えば、指数関数であればよい。指数関数が用いられる場合には、そのパラメータである底が減衰速度の指標として得られる。回帰分析の手法として、線形分析が用いられてもよいが、これに限らず、非線形分析が用いられてもよい。
 以上のように、本実施形態に係る生体組織の評価方法は、生体組織のうち、特定の一部分の画像を短時間に多数回取得し、スペックルバリアンス、OCT相関の減衰速度等を算出することによって、生体組織の活動性を定量評価するものである。これにより、従来は見過ごされていた微小な変化(ゆらぎ)を検出することができ、生体組織のダイナミクスの定量的な評価を実現することができる。
 以下、実施例により本実施形態の効果をより明らかなものとする。なお、本実施形態は、以下の実施例に限定されるものではなく、その要旨を変更しない範囲で適宜変更して実施することができる。
(比較例1)
 人がん組織由来の細胞を球状に培養した腫瘍凝集体サンプルに対して、4時間おきに、従来のOCT撮影を行った。
(実施例1)
 比較例1と同様の腫瘍凝集体サンプルに対して、一例として、2時間おきに、13ms間隔で100回の高速OCT撮影を行った。観測期間間隔を2時間としたのは、試料とする腫瘍凝集体組織が生体から取り出してから死に至るまでの期間(典型的には、1~3日程度)における活動状態の変化傾向を捉えるために十分なためである。また、1回の観測期間を13ms間隔で100回、即ち、1.3秒としたのは、この観測期間が観測期間間隔よりも十分に短く、組織を構成する細胞の運動もしくは細胞内の活動(例えば、数Hz~20Hz)に応じた光学的特性の時間変化を把握するのに十分なためである。続いて、各回のOCTの画像から得られるOCT信号強度の時間変化に基づいて、スペックルバリアンスを算出した。
(実施例2)
 比較例1と同様の腫瘍凝集体サンプルに対して、実施例1と同様の高速OCT撮影を行った。続いて、各回のOCTの画像から得られるOCT信号強度の時間変化に基づいて、OCT相関の減衰速度を算出した。
(実施例3)
 比較例1と同様の腫瘍凝集体サンプルに対して、実施例1と同様の高速OCT撮影を行った。続いて、各回のOCTの画像から得られるOCT信号強度に基づいて、深さ方向に対するOCT信号強度の減衰係数(AC:Attenuation Coefficient)を算出した。
 図2は、比較例1、実施例1、2の方法を適用して得られた腫瘍凝集体サンプルの画像である。上段、中段、下段の画像が、それぞれ比較例1、実施例1、2に対応している。0hr、8hr、28hrとの数値は、それぞれ人体から切り出され培養を開始した時点からの経過時間を示す。比較例1の画像では、腫瘍凝集体サンプルの状態に、時間経過による変化は見られない。これに対し、実施例1、2の画像では、初期の段階で中心に暗い部分があり、その周囲に明るい部分あるが、時間経過により、明るい部分が暗く変わってゆく様子が見られる。暗い部分が死んでいる状態を表し、明るい部分が生きている状態を表していると考えられるため、これらの画像の変化から、腫瘍凝集体サンプルの生死状態について確認することができる。
 図3A、図3Bは、それぞれ図2の実施例1、2の画像に示すサンプルから個々の観測期間について算出されたスペックルバリアンス、OCT相関の減衰速度の経時変化を示す図である。図3Cは、実施例3の方法を適用して得られるOCT信号強度の減衰係数を示す図である。図3A、図3B、図3Cの横軸は、いずれも腫瘍凝集体サンプルを切り出してからの経過時間を示しており、縦軸は、それぞれ、腫瘍凝集体サンプル全体でのスペックルバリアンス、OCT相関の減衰速度、OCT信号強度の減衰係数の平均値を示している。いずれの図においても、時間経過とともに、それぞれの平均値の減少傾向が見られる。これら3つの図の比較から、OCT相関の減衰速度の傾きが大きいため、腫瘍サンプル定量評価する上で最も適していると考えられる。そして、スペックルバリアンスに対する傾きの大きさが次に大きく、従来提案されたOCT信号強度の減衰係数に対する傾きよりも有意に大きい。なお、A.U.とは、任意単位(arbitrary unit)を示す。
 なお、複数回の前記OCTの画像のうち1回目のものは、それ以降のサンプルとの相関係数の減衰傾向を正しく反映せず、その減衰傾向から推定される推定値と有意差を生ずるため、定量評価においては、時間シフト量τがゼロとなる場合における相関係数を除外し、ゼロ以外となる非零の時間シフト量τに対する相関係数を用いてOCT相関の減衰速度を算出することが好ましい。
 図4A、図4Bは、それぞれ実施例1、2の方法を適用した生体組織の生死状態を示す図である。図4A、図4Bの横軸は、いずれも経過時間を示しており、縦軸は、生きている細胞または死んでいる細胞の含有率を示している。生死の境界について、実施例1(SV)では3.0とし、実施例2(OCDS)では5×10-4ms-1とした。いずれのグラフにおいても、時間経過とともに、生きている状態の細胞が減少し、死んでいる細胞が増加する傾向が見られる。
(比較例2)
 生体組織の一部として切り出した、ネズミの肝臓サンプルに対して、1時間おきに、従来のOCT撮影を行った。
(実施例4)
 比較例2と同様の肝臓サンプルに対して、1時間おきに、10ms間隔で100回の高速OCT撮影を行った。続いて、各回のOCTの画像から得られるOCT信号強度の時間変化に基づいて、スペックルバリアンスを算出した。
(実施例5)
 比較例2と同様の肝臓サンプルに対して、実施例4と同様の高速OCT撮影を行った。続いて、各回のOCTの画像から得られるOCT信号強度の時間変化に基づいて、OCT相関の減衰速度を算出した。
(実施例6)
 比較例2と同様の肝臓サンプルに対して、実施例4と同様の高速OCT撮影を行った。続いて、各回のOCTの画像から得られるOCT信号強度に基づいて、深さ方向に対するOCT信号強度の減衰係数(AC)を算出した。
 図5は、比較例2、実施例4、5の方法を適用して得られた肝臓サンプルの画像である。上段、中段、下段の画像が、それぞれ比較例2、実施例4、5に対応している。0hr、8hr、16hrなどの数値は、それぞれ母体とするネズミから切り出された時点からの経過時間を示す。比較例2の画像では、肝臓サンプルの状態に、時間経過による変化は見られない。これに対し、実施例4、5の画像では、初期の段階で下側に暗い部分があり、上側に明るい部分あるが、時間経過により、明るい部分が暗く変わってゆく様子が見られる。暗い部分が死んでいる状態を表し、明るい部分が生きている状態を表していると考えられるため、これらの画像の変化から、肝臓サンプルの生死状態について確認することができる。
 図6A、図6Bは、それぞれ実施例4、5の画像に示すサンプルから個々の観測期間について算出されたスペックルバリアンス、OCT相関の減衰速度を示す図である。図6Cは、実施例6の方法を適用して得られたOCT信号強度の減衰係数を示す図である。図6A、図6B、図6Cの横軸は、いずれも肝臓サンプルを切り出してからの経過時間を示しており、縦軸は、それぞれ、肝臓サンプル全体でのスペックルバリアンス、OCT相関の減衰速度、減衰係数の平均値を示している。いずれの図においても、時間経過とともに、それぞれの平均値について、傾きが異なる二段階の減少傾向が見られる。これら3つの図の比較から、OCT相関の減衰速度の傾きが大きいため、腫瘍サンプル定量評価する上で最も適していると考えられる。そして、スペックルバリアンスに対する傾きの大きさが次に大きく、従来提案されたOCT信号強度の減衰係数に対する傾きよりも有意に大きい。
 なお、複数回の前記OCTの画像のうち1回目のものは、それ以降のサンプルとの相関係数の減衰傾向を正しく反映せず、その減衰傾向から推定される推定値と有意差を生ずるため、定量評価においては、τがゼロとなる場合を除外して、OCT相関の減衰速度を算出することが好ましい。
 図7A、図7Bは、それぞれ実施例4、5の方法を適用した生体組織の生死状態を示す図である。図7A、図7Bの横軸は、いずれも経過時間を示しており、縦軸は、生きている細胞または死んでいる細胞の含有率を示している。生死の境界について、実施例4(SV)では3.0とし、実施例5(OCDS)では5×10-4ms-1とした。いずれのグラフにおいても、時間経過とともに、生きている状態の細胞が減少し、死んでいる細胞が増加する傾向が見られる。このように、本実施形態によれば各観測期間においてスペックルバリアンス(SV)、OCT相関(OCDS)の減衰速度などの時間変動特性値を算出して、これらをより長期間にわたり観測をすることで細胞の生死に関するダイナミクスを観察することができる。
<第2の実施形態>
 次に、本発明の第2の実施形態について説明する。
 上記のダイナミックOCTと呼ばれる解析手法によれば、細胞など生体組織をなす内因性散乱因子の局所的活動を表現可能とする。しかしながら、非特許文献1に記載の手法では、表現される局所的活動の定量性に乏しいため、生体の活性を正しく評価することが困難なことがあった。
 他方、糖尿病研究、循環器研究などの応用面では、微小血管、リンパ管、心筋など特定の組織に特異性を有する可視化するための手法としてダイナミクスイメージングの応用が期待されている。生体組織の活性を定量的に評価するために、生体組織の同じ位置において、複数回のOCT計測を行い、各回の計測により得られるOCT信号強度の時間変化を解析することも考えられる。しかしながら、組織内のあらゆる微小な動きがOCT信号の時間変化に反映されるため、単に信号強度の時間変化を解析するだけでは、特定の組織に対する活性を定量的に評価することができないことがあった。本実施形態は、この点に鑑みて提案されたものである。
 図9は、本実施形態に係るOCTシステム1の一例を示す構成図である。OCTシステム1は、PS-OCTを構成する。PS-OCTは、試料Smに対して偏光状態が既知である入射光を照射し、試料Smから反射した反射光と参照光とを干渉させた干渉光を取得するための光学系を備える。また、OCTシステム1は、光学系により取得された干渉光の偏光状態から、試料Smにおける偏光状態の変化特性を解析する計測信号処理装置200を備える。計測信号処理装置200は、OCT信号を用いて試料Smとする生体組織の状態を解析する評価装置として機能する。計測信号処理装置200は、解析した変化特性を可視化した画像を生成する。
 試料Smとする観測対象の物体は、主に、人間もしくは動物などの生体の一部である。より具体的には、眼底、血管、歯牙、皮下組織などのいずれであってもよい。これにより、試料Sm内部の状態を非侵襲で計測または観測することができる。そのため、生体内組織、例えば、眼底などの微小血管、リンパ管、心筋などの診断への応用が期待されている。
 図9に例示されるOCTシステム1は、光源102で発生する光の波長を掃引してスペクトル干渉信号を得るための波長掃引型OCT(SS-OCT:Swept Source-OCT)を応用した観測システムをなす。OCTシステム1は、光源102から出射される光をプローブアーム(後述)と参照アーム130に分岐して入射させる。OCTシステム1は、プローブアームに分岐した光を水平偏光成分と、垂直偏光成分に分離し、相互間で光路長が異なる偏光成分を含む光を計測対象の試料Smに走査(B-スキャン)しながら照射し、試料Smから反射した反射光(物体光)を取得する。OCTシステム1は、参照アーム130に分岐した参照光と試料Smから反射、散乱または、その両者により得られた成分である反射光を干渉させて干渉光を取得する。なお、本願では、試料Smに対して光が照射される方向を深さ方向とする。試料Smの深さ方向への観測点の走査による計測信号の取得は、A-スキャンと呼ばれる。SS-OCTでは、波長掃引光源を用いることでA-スキャンが実現される。B-スキャンとは、試料Smの深さ方向に垂直な方向への走査を指す。
 OCTシステム1は、光源102、カプラー104、偏光遅延ユニット110、サーキュレータ120、プローブ128、参照アーム130、偏光分離検出ユニット150、光検出器190および計測信号処理装置200を備える。光源102、カプラー104、偏光遅延ユニット110、サーキュレータ120、プローブ128、偏光分離検出ユニット150および光検出器190は、それぞれ光学系を構成する構成要素である。光学系の構成要素間は光路として光ファイバを用いて結合されている。
 光源102は、周期的に所定の波長幅(例えば、40~120nm)内で掃引する波長を有する光を発生させる波長掃引光源(Wavelength Swept Source)である。光源102は、例えば、SLD(スーパールミネセントダイオード;Superluminescent Diode)など近赤外の波長(例えば、1000~1400nm)を有する。光源102が発生させた光は、カプラー104に入射する。
 カプラー104は、光源102から入射される光を所定の強度比でプローブアームと参照アーム130の2系統に分離する。プローブアームへの光強度と参照アーム130への光強度の比は、例えば、90%:10%である。プローブアームに供給された光は偏光分離検出ユニット(PDDU:Polarization Diversity Detection Unit)150に供給される。プローブアームは、ファイバコリメータ106、偏光制御器108、偏光遅延ユニット(PDU:Polarization Delay Unit)110、サーキュレータ120、ファイバコリメータ122、偏光制御器124、対物レンズ126およびプローブ128がその順序で接続されてなる経路である。プローブアームは、サンプルアームまたは測定アームとも呼ばれる。プローブアームに供給される光は、ファイバコリメータ106、偏光制御器108を経由して偏光遅延ユニット110に入射される。他方の系統の光は、参照アーム130を経由してPPDU150に入射される。なお、偏光制御器108は、入射光の強度を所定の十分な強度に増幅し、増幅した光を出射する。
 PDU110は、直線偏光器(Linear Polarizer)112、偏光ビームスプリッタ(PBS:Polarization Beam Splitter)114および2個の直角プリズム(RAP:Right Angle Prism)116、118を備える。PDU110は、入射光を相互に直交する2つの偏光状態を有する成分として水平偏光成分と垂直偏光成分に分離し、分離された各成分を合波して得られる光をサーキュレータ120に供給する。
 直線偏光器112は、カプラー104から入射される光の偏光状態を直線偏光に変換し、変換した光をPBS114に出射する。出射される光の水平偏光成分と垂直偏光成分を等しくするため、直線偏光器112の偏光角度を45°と設定されている。PBS114は、その表面が入射角を45°とする方向に配置された反射層を有し、直線偏光器112から反射層に入射される入射光のうち垂直偏光成分が透過光として透過し、反射層の表面を水平偏光成分が反射光として反射する。PBS114からの水平偏光成分を含む反射光、垂直偏光成分を含む透過光は、それぞれRAP116、118に入射される。他方、PBS114の反射層は、RAP116から反射層に入射される水平偏光成分を含む入射光を透過した透過光と、PBS118から入射される垂直偏光成分を含む入射光に対する反射光を合波し、合波された光をサーキュレータ120に出射する。
 RAP116、118は、それぞれ光路に平行な断面が直角二等辺三角形となる形状を有し、直角二等辺三角形の底辺がPBS114からの光路に直交する向きに配置される。PBS114から入射される光は、直角二等辺三角形の底辺に平行な側面を透過し、底辺に向かい合う2つの辺のうちの一辺に平行な側面で反射され、その反射光を他方の辺に平行な側面で反射されて底辺に平行な側面に戻り、その側面を透過してPBS114に入射される。なお、PBS114とRAP116との光路長と、PBS114とRAP118との光路長が有意に異なるようにRAP118の位置を予め調整しておく。これにより、PBS114から試料Smに入射される水平偏光成分と垂直偏光成分が、相互間で所定の位相差をもって重畳して出射される。
 サーキュレータ120は、PDU110から入射される光をファイバコリメータ122および偏光制御器124を経由して対物レンズ126に出射する。対物レンズ126は、自部に入射される光を集光し、プローブ128を経由して試料Smに照射する。試料Smから反射、散乱、または、その両者により取得された光はプローブ128を経由して対物レンズ126で平行光に変換され計測ビームとして偏光制御器124およびファイバコリメータ122を経由してサーキュレータ120に戻り、PPDU150に入射される。
 参照アーム130は、ファイバコリメータ132、ファイバブラッググレーティング(FBG:Fiber Bragg Grating)134、ファイバコリメータ136、ディレイライン(Delay Line)138、ファイバコリメータ140および偏光制御器142がその順序で接続されてなる経路である。
 FBG134は、入射される光のうち特定の波長の成分を反射光として反射し、それ以外の残りの成分を透過して透過光としてファイバコリメータ136を経由してディレイライン138に入射される。FBG134からの反射光は、ファイバコリメータ132を経由してカプラー104に戻り、カプラー104から光検出器190に入射される。FBG134から反射される成分の帯域は、光源102で生成される光の帯域よりも十分に狭い。光検出器190は、FBG136からの反射光の強度を検出し、検出した強度を示す強度信号をトリガ信号として計測信号処理装置200に出力する。トリガ信号は、A-スキャンのトリガとして用いられる。光源102で生成される光の波長が所定の波長幅の範囲内で周期的に変化するが、波長が所定の波長となるタイミングが光検出器190で検出され、そのタイミングでA-スキャンが光学系制御部212でリセットされる。光検出器190には、例えば、その波長幅の下限を検出対象の波長として予め設定しておく。SS-OCTでは、観測目標とする観測点の深度がプローブ光の波長により定まるためである。
 ディレイライン138は、FBG134から入射される入射光を遅延させ、遅延させた光をファイバコリメータ140および偏光制御器142を経由してPPDU150に出射する。ディレイライン138は、入射光に対する遅延量を可変とし、プローブアームの光路長と参照アーム130の光路長が等しくなるように遅延量を調整させておく。なお、偏光制御器142は、入射光の強度を所定の強度に調整し、強度が調整された光を出射する。
 PPDU150は、直線偏光器152、非偏光ビームスプリッタ(NPBS:Non-Polarization Beam Splitter)154、2個のPBS156、158、4個の受光器162、164、166、168および2個のバランス偏光検出器(BPD:Balanced Polarization Detector)170、172を備える。
 直線偏光器152は、参照アーム130から入射される光の偏光状態を直線偏光に変換し、変換した光をNPBS154に出射する。直線偏光器152の偏光角度は45°と設定されている。NPBS154は、参照アーム130から直線偏光器152を経由して入射される入射光と、プローブアームから入射される入射光を合波する。NPBS154は、参照アーム130からの入射光、プローブアームからの入射光のそれぞれに対して、その表面が入射角を45°とする方向に配置された反射層を有する。反射層は、参照アーム130からの入射光を透過して得られる透過光と、プローブアームからの入射光を反射して得られる反射光を合波し、合波により得られる干渉光をPBS158に出射する。反射層は、参照アーム130からの入射光を透過して得られる反射光と、プローブアームからの入射光を透過して得られる透過光を合波し、合波により得られる干渉光をPBS156に出射する。
 PBS156は、NPBS154から入射される干渉光を水平偏波成分と垂直偏波成分に分離し、分離した水平偏波成分と垂直偏波成分をそれぞれ受光器162、166に出射する。受光器162、166は、それぞれPBS156から入射される水平偏波成分、垂直偏波成分を受光し、それぞれ第1水平偏波成分、第1垂直偏波成分としてBPD170、172に導光する。第1水平偏波成分、第1垂直偏波成分は、それぞれ試料Smへ入射される水平偏波成分に基づく干渉光の水平偏波成分、垂直偏波成分に相当する。
 PBS158は、NPBS154から入射される光を水平偏波成分と垂直偏波成分に分離し、分離した水平偏波成分と垂直偏波成分をそれぞれ受光器164、168に出射する。受光器164、168は、それぞれPBS158から入射される水平偏波成分、垂直偏波成分を受光し、それぞれ第2水平偏波成分、第2垂直偏波成分としてBPD170、172に導光する。第2水平偏波成分、第2垂直偏波成分は、それぞれ試料Smへ入射される垂直偏波成分に基づく干渉光の水平偏波成分、垂直偏波成分に相当する。
 BPD170は、受光器162、166からそれぞれ導光される第1水平偏波成分、第2水平偏波成分を検出し、検出した第1水平偏波成分、第2水平偏波成分の強度を示すアナログの電気信号である第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号に変換する。BPD170は、生成した第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号を低域通過フィルタ(LPF:Low Pass Filter)182と高域通過フィルタ(HPF:High Pass Filter)186を経由して計測信号処理装置200に出力する。
 BPD172は、受光器164、168からそれぞれ導光される第1垂直偏波信号、第2垂直偏波成分を検出し、検出した第1垂直偏波成分、第2垂直偏波成分の強度を示すアナログの電気信号である第1垂直偏波スペクトル干渉信号、第2垂直偏波スペクトル干渉信号に変換する。BPD172は、生成した第1垂直偏波スペクトル干渉信号、第2垂直偏波スペクトル干渉信号をLPF184とHPF188を経由して計測信号処理装置200に出力する。第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号、第1垂直偏波スペクトル干渉信号、第2垂直偏波スペクトル干渉信号は、それぞれ各1フレームのOCT画像の生成に用いられ、ジョーンズ行列断層撮影(JM-OCT:Jones Matrix-OCT)を実現可能とする。
 なお、光源102が発光する光の波長の幅をBPD170、172は、1回のA-スキャンごとに所定のサンプル数(例えば、400~2000サンプル)の信号値を所定のサンプリング周波数でサンプリングする。LPF182、184、HPF186、188は、例えば、それぞれチェビシェフフィルタである。LPF182、184のカットオフ周波数は、例えば、62MHzである。HPF186、188のカットオフ周波数は、例えば、1MHzである。
(計測信号処理装置)
 次に、本実施形態に係る計測信号処理装置200の構成例について説明する。図10は、本実施形態に係る計測信号処理装置200の構成例を示すブロック図である。計測信号処理装置200は、PPDU150から入力される第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号、第1垂直偏波スペクトル干渉信号および第2垂直偏波スペクトル干渉信号から試料内の観測点における偏光特性をそれぞれ解析し、解析した偏光特性を示す偏光特性値の時間変更特性を解析し、解析した時間変動特性を示す時間変動特性値を定める。計測信号処理装置200は、観測点ごとに定めた時間変動特性値に応じた色または階調を示す信号値に変換し、変換した信号値を観測点に対応する画素ごとに画像データを生成し、生成した画像データを出力してもよい。
 計測信号処理装置200は、制御部210と、記憶部230と、入出力部240と、表示部250と、操作部260と、を含んで構成される。制御部210の一部または全部の機能は、例えば、CPU(Central Processing Unit)等のプロセッサを含んで構成されるコンピュータとして実現される。プロセッサは、予め記憶部230に記憶させておいたプログラムを読み出し、読み出したプログラムに記述された指令で指示される処理を行って、その機能を奏する。本願では、プログラムに記述された指令で指示される処理を行うことを、プログラムを実行する、プログラムの実行、などと呼ぶことがある。制御部210の一部または全部は、プロセッサなどの汎用のハードウェアに限られず、LSI(Large Scale Integration)、ASIC(Application Specific Integrated Circuit)等の専用のハードウェアを含んで構成されてもよい。
 制御部210は、光学系制御部212、計測信号取得部214、偏光解析部216、変動特性解析部218、画像処理部220および出力処理部222を含んで構成される。本実施形態に係る計測信号取得部214と偏光解析部216は、OCT信号を取得し、取得したOCT信号に基づく信号値を試料の観測点ごとに取得する計測部として機能する。変動特性解析部218は、所定期間ごとに、取得された信号値の時間変動特性を示す時間変動特性値を算出する評価部として機能する。
 光学系制御部212は、プローブの位置を可変にする駆動機構を駆動させ、試料Smの観測点を走査する(B-スキャン)。光学系制御部212は、試料Smの観測点を試料Smの深さ方向(以下、z方向)に交差する所定の方向(例えば、試料Smの深さ方向に直交する正面上のx方向)に走査される。光学系制御部212は、例えば、光検出器190から入力されるトリガ信号が所定の強度以上となるとき、その方向へのプローブ128の位置の移動を示す制御信号を駆動機構に出力する。この移動距離は、予め設定されたx方向の観測点間隔に相当する。また、このとき、光学系制御部212は、基準点からのx方向の観測点数が所定の個数(ライン数)に達するごとに、プローブ128の位置を基準点に戻す。その後、計測信号取得部214は、次のフレームの計測信号を取得する。計測信号の繰り返しにより、時間経過に応じて計測信号がフレームごとに累積される。図12に示す例では、x方向、z方向、時刻tが、それぞれ右方、下方、右上方に示され、各フレームが個々の長方形で示される。
 計測信号取得部214には、PPDU150から入出力部240を経由して第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号、第1水平偏波スペクトル干渉信号および第2垂直偏波スペクトル干渉信号が計測信号として入力される。計測信号取得部214は、第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号、第1垂直偏波スペクトル干渉信号および第2垂直偏波スペクトル干渉信号をそれぞれフーリエ変換し、それぞれの複素振幅を示す第1水平偏波OCT信号、第2水平偏波OCT信号、第1垂直偏波OCT信号および第2垂直偏波OCT信号を観測点ごとに算出する。計測信号取得部214は、算出した第1水平偏波OCT信号、第2水平偏波OCT信号、第1垂直偏波OCT信号および第2垂直偏波OCT信号を偏光解析部216に出力する。
 偏光解析部216は、計測信号取得部214から入力される第1水平偏波OCT信号、第2水平偏波OCT信号、第1垂直偏波OCT信号および第2垂直偏波OCT信号に基づいて観測点における偏光特性を示す偏光特性値を試料Sm内の観測点ごとに算出する。偏光解析部216は、計測に係る光波の偏光特性を示すジョーンズ行列を構成し、構成したジョーンズ行列から偏光特性を表す指標値として所定の偏光特性値を算出する。
 ジョーンズ行列は、偏光状態の変化を示す2行2列の行列である。ジョーンズ行列の第1列、第2列には、それぞれ第1ジョーンズベクトル、第2ジョーンズベクトルが配置される。第1ジョーンズベクトル、第2ジョーンズベクトルは、相互に直交する偏光成分を有する入射光を用いて取得される。即ち、第1ジョーンズベクトルは、第1水平偏波OCT信号と第1垂直偏波OCT信号でそれぞれ示される複素振幅を要素として含む2次元のベクトルである。第2ジョーンズベクトルは、第2水平偏波OCT信号と第2垂直偏波OCT信号でそれぞれ示される複素振幅を要素として含む2次元のベクトルである。ジョーンズ行列の各列には入射に係る偏光状態、各行には検出された偏光状態が対応する。偏光解析部216は、算出した観測点ごとの偏光特性値を変動特性解析部218に出力する。なお、本願では、この段階でOCT信号またはその周波数領域の変換係数から直接構成されるジョーンズ行列を計測ジョーンズ行列(Measured Jones Matrix)と呼ぶことがある。
 変動特性解析部218は、偏光解析部216から入力される偏光特性値について予め設定された所定期間である観測期間における時間変化特性を示す指標値として所定の時間変動特性値を算出する。観測期間は、例えば、フレームレートが60フレーム/秒であるとき、典型的には150-600フレーム程度である。変動特性解析部218は、算出した時間変動特性値を画像処理部220に出力する。なお、偏光特性値と時間変動特性値の例については、後述する。
 画像処理部220は、変動特性解析部218から入力される時間変動特性値を、所定の関数を用いて画素ごとのビット深度で表現可能な所定の値域内の画素値に変換する。画像処理部220は、例えば、時間変動特性値に対するシグモイド関数(Sigmoid Function)の関数値に所定の倍数を乗じて得られる乗算値に所定のオフセット値を加算して画素値を算出する。時間変動特性値を画素値に変換するための関数は、シグモイド関数に限られず、一次関数、対数関数など、入力値に対して得られる関数値が、入力値の増加に対して関数値が単調に増加または減少する関数であれば利用可能である。画像処理部220は、観測点ごとに変換した画素値を示す出力画像データを生成し、生成した出力画像データを表示部250に出力する。
 画像処理部220は、出力処理部222から入力される制御信号に応じて出力画像データを記憶部230に記憶してもよい。
 出力処理部222は、操作部260から入力される操作信号に基づいて表示画像を示す出力画像データの生成または出力を制御する。操作信号により、例えば、表示画像の表示または記憶の要否、観測対象領域などがパラメータとして指示される。制御部210が、複数の種別の時間変動特性値を算出する能力を有する場合には、それらの複数の種別の時間変動特性値のうち、算出もしくは表示対象の種別が操作信号により指示されてもよい。出力処理部222は、操作により設定可能とするパラメータ、パラメータの設定、設定可能とするパラメータを案内するための設定画面を表示部に表示させ、画像の表示に係るユーザインタフェースを構成してもよい。
 例えば、表示画像の表示の要否を示す操作信号が入力されるとき、出力処理部222は、その表示の要否を示す制御信号を画像処理部220に出力する。画像処理部220は、出力処理部222から表示要を示す制御信号が入力されるとき出力画像データを表示部に出力し、出力処理部222から表示否を示す制御信号が入力されるとき出力画像データを表示部に出力しない。
 観測対象領域を示す操作信号が入力されるとき、出力処理部222は、その観測対象領域内のx座標またはy座標の範囲を示す制御信号を計測信号取得部214に出力する。計測信号取得部214は、出力処理部222から入力される制御信号で示される範囲でB-スキャンを実行する。観測対象領域は、例えば、試料Smの表面の範囲がパラメータとして定められる。
 記憶部230は、上記のプログラムの他、制御部210が実行する処理に用いられる各種のデータ、制御部210が取得した各種のデータを記憶する。記憶部230は、例えば、ROM(Read Only Memory)、フラッシュメモリ、HDD(Hard Disk Drive)などの不揮発性の(非一時的)記憶媒体を含んで構成される。記憶部230は、例えば、RAM(Random Access Memory)、レジスタなどの揮発性の記憶媒体を含んで構成される。
 入出力部240は、無線または有線で他の機器と各種のデータを入出力可能に接続する。入出力部240は、例えば、入出力インタフェースを備える。入出力部240は、例えば、偏光分離検出ユニット150と光検出器190に接続される。
 表示部250は、制御部210から入力される出力画像データに基づく画像を表示する。表示部250は、例えば、液晶ディスプレイ、有機エレクトロルミネッセンスディスプレイなどのいずれを備えてもよい。
 操作部260は、例えば、ボタン、つまみ、ダイヤル、マウス、ジョイスティックなど、ユーザの操作を受け付け、受け付けた操作に応じた操作信号を生成する部材を含んで構成されてもよい。操作部260は、取得した操作信号を制御部210に出力する。操作部260は、他の機器(例えば、リモートコントローラ等の携帯機器)から操作信号を無線または有線で受信する入力インタフェースであってもよい。
(OCT信号処理)
 次に、本実施形態に係るOCT信号処理の例について説明する。図11は、本実施形態に係るOCT信号処理の一例を示すフローチャートである。
(ステップS102)計測信号取得部214は、光学系から第1水平偏波スペクトル干渉信号、第2水平偏波スペクトル干渉信号、第1垂直偏波スペクトル干渉信号および第2垂直偏波スペクトル干渉信号を取得し、取得された信号のそれぞれから第1水平偏波OCT信号、第2水平偏波OCT信号、第1垂直偏波OCT信号および第2垂直偏波OCT信号を算出する。
(ステップS104)偏光解析部216は、計測信号取得部214により取得された第1水平偏波OCT信号、第2水平偏波OCT信号、第1垂直偏波OCT信号および第2垂直偏波OCT信号に基づいて観測点ごとにジョーンズ行列を構成する。偏光解析部216は、構成したジョーンズ行列から所定の偏光特性値を算出する。
(ステップS106)変動特性解析部218は、予め設定された観測期間における偏光解析部216により算出された変更特性値から所定の時間変動特性値を算出する。
(ステップS108)画像処理部220は、変動特性解析部218で算出された観測点ごとの時間変動特性値を、個々の観測点に対応する画素における画素値に変換する。
(ステップS110)画像処理部220は、変換した画素値を示す出力画像データを表示部250に出力する(画像表示)。よって、表示部250は、観測対象領域における偏光特性を示す時間変動特性値の観測対象領域における分布を可視化する。
(偏光特性値、時間変動特性値の例)
 次に、偏光特性値と時間変動特性値の例について説明する。偏光解析部216は、偏光特性値として、例えば、位相遅延(Phase Retardation)を算出する。位相遅延は、複屈折(Birefringence)により生ずる通常光線と異常光線との位相差である。即ち、通常光線と異常光線は、試料の光学軸に対して互いに直交する偏光方向を有し、それぞれ異なる速度で試料を透過する。位相遅延のうち表面から試料内の注目する観測点まで累積された位相遅延は累積位相遅延(CPR:Cumulative Phase Retardation)と呼ばれる。偏光解析部216は、(3)式に示すように深さzの観測点に対応する試料表面上の位置におけるジョーンズ行列の逆行列Jm0 -1に深さzの観測点における計測ジョーンズ行列Jmzを乗じて観測点における累積ジョーンズ行列(Cumulative Jones Matrix)Jczを算出する。累積ジョーンズ行列Jczは試料表面から深さzの観測点までの偏光状態の変化を示す。
Figure JPOXMLDOC01-appb-M000003
 偏光解析部216は、累積ジョーンズ行列Jczに対して固有値分解を行って得られる2個の固有値λc1、λc2の位相差arg(λc1λc2 )をCPRとして算出することができる。なお、arg(…)は、複素数…の偏角を示し、~は、複素数~の複素共役を示す。
 なお、偏光解析部216は、偏光特性値の他の例として、局所的な深さ領域における偏光位相遅延である局所位相遅延(LPR:Local Phase Retardation)を算出してもよい。より具体的には、偏光解析部216は、(4)式に示すように深さzの観測点における累積ジョーンズ行列の逆行列Jcz1 -1に深さzの観測点における累積ジョーンズ行列を、Jcz2を乗じて得られる深さz、zの観測点間における局所ジョーンズ行列Jl12を算出する。局所ジョーンズ行列Jl12は深さzの観測点から深さzの観測点までの偏光状態の変化を示す。
Figure JPOXMLDOC01-appb-M000004
 深さzの観測点から深さzの観測点までの厚みΔz12は、例えば、操作部260から操作信号で指示される観測対象の組織の厚みであってもよいし、z方向の観測点間距離であってもよい。偏光解析部216は、局所ジョーンズ行列Jl12の2個の固有値λl1、λl2の位相差arg(λl1λl2 )をLPRとして算出することができる。なお、偏光解析部216は、深さzの観測点における計測ジョーンズ行列の逆行列Jmz1 -1に深さzの観測点における計測ジョーンズ行列Jmz2を乗じて深さz、zの観測点間における局所ジョーンズ行列Jl12を算出してもよい。
 また、偏光解析部216は、偏光特性値の他の例として、(5)式に示すように、LPRを2kΔz12で除算して深さz、zの観測点間における複屈折率b12を算出してもよい。但し、kは、入射光の中心波長を示す。
Figure JPOXMLDOC01-appb-M000005
 変動特性解析部218は、所定の観測期間内における時間変動特性値の例として、偏光特性値の分散を算出する。算出される時間変動特性値は、偏光特性値で示される偏光状態の時間変動の大きさを示す。(6)式において、σ (x,z)はLPRの分散を示す。N、φ(x,z,t)、<φ(x,z)>は、それぞれ所定の観測期間におけるフレーム数、LPR、LPRの時間平均値を示す。(x,z)は観測点のx-z平面内の座標を示す。tは、第i番目のサンプル時刻を示す。
Figure JPOXMLDOC01-appb-M000006
 式(6)は、LPRの分散を例示するが、変動特性解析部218は、LPRに代えて、CPRの分散または複屈折率の分散を算出してもよい。変動特性解析部218は、これらの分散の平方根を標準偏差として算出してもよい。
 変動特性解析部218は、時間変動特性値の他の例として、偏光特性値の対数値の分散または標準偏差を算出してもよい。分散または標準偏差を算出する過程で、各時刻の偏光特性値の対数値から時間平均値の対数値が差し引かれるため、偏光特性値に潜在的に乗じられた定数が消去される。従って、実質的な時間変動特性の評価に適用されうる。また、対数値をとることで、スケールが大きく異なる現象同士の比較に適用されうる。(7)式において、logσ (x,z)はLPRの対数値の分散を示す。
Figure JPOXMLDOC01-appb-M000007
 (7)式は、LPRの対数値の分散を例にするが、変動特性解析部218は、LPRの対数値の分散に代えて、CPRの対数値の分散または複屈折率の対数値の分散を算出してもよい。変動特性解析部218は、これらの分散の平方根を標準偏差として算出してもよい。
 変動特性解析部218は、所定の観測期間内における時間変動特性値として、例えば、偏光特性値の動的コントラスト(Dynamic Contrast)を算出してもよい。偏光特性値の動的コントラストは、偏光特性値に対する標準偏差を時間平均値で除算することによって正規化された値に相当する。(8)式において、φはLPRの動的コントラストを示す。
Figure JPOXMLDOC01-appb-M000008
 (8)式は、LPRの動的コントラストを例にするが、変動特性解析部218は、LPRに代えて、CPRまたは複屈折率の動的コントラストを算出してもよい。動的コントラストを算出する過程では、標準偏差が時間平均値で正規化されるため、動的コントラストは標準偏差よりも実質的な時間変動特性の評価に適用されうる。
 なお、変動特性解析部218は、時間変動特性値の他の例として、時間偏光均一度(TPU:Temporal Polarization Uniformity)を算出してもよい。図18を参照しながら、TPUの算出方法について説明する。まず、偏光解析部216は、偏光特性値の他の例として、観測点に係る計測ジョーンズ行列Jmzの部分空間をなす第1ジョーンズベクトルJ、第2ジョーンズベクトルJを、それぞれ第1ストークスベクトルS、第2ストークスベクトルSに変換する(図18、ステップS122)。(9)式に示すように第1ジョーンズベクトルJ、第2ジョーンズベクトルJは、それぞれジョーンズ行列Jmzの第1列の要素、第2列の要素を含んで構成される。
Figure JPOXMLDOC01-appb-M000009
 第1ストークスベクトルS、第2ストークスベクトルSは、式(10)に示すように、それぞれ第1ジョーンズベクトルJ、第2ジョーンズベクトルJで示される偏光状態を表現する4次元のベクトルであり、それぞれ第1ジョーンズベクトルJ、第2ジョーンズベクトルJの要素値を含んで構成される。以下の説明では、第1ストークスベクトルS、第2ストークスベクトルSのそれぞれをなす各4個の要素値s10~s13、s20~s23を第0~第3ストークスパラメータs10~s13、s20~s23と呼ぶ。
Figure JPOXMLDOC01-appb-M000010
 (10)式に示すように、第0ストークスパラメータs10、s20は、水平方向成分のパワー|g1H、|g2Hと垂直方向成分のパワー|g1V、|g2Vの和、つまり光の全体の強度を示す。第1ストークスパラメータは、水平方向成分のパワー|g1H、|g2Hから垂直方向成分のパワー|g1V、|g2Vの差、つまり、相互に直交する成分間の差を示す。第2ストークスパラメータは、水平方向成分と垂直方向成分の複素共役の積g1H1V 、g2H2V の実部の2倍、つまり、水平方向成分g1H、g2Hと垂直方向成分g1V、g2Vそれぞれの強度の積|g1H||g1V|、|g2H||g2V|、水平方向成分、と垂直方向成分の位相差δ、δの余弦値cosδ、cosδに2を乗じて得られる値に相当する。第3ストークスパラメータは、水平方向成分と垂直方向成分の複素共役の積g1H1V 、g2H2V の虚部の2倍、つまり、水平方向成分g1H、g2Hと垂直方向成分g1V、g2Vそれぞれの強度の積|g1H||g1V|、|g2H||g2V|、水平方向成分、と垂直方向成分の位相差δ、δの正弦値sinδ、sinδに2を乗じて得られる値に相当する。
 変動特性解析部218は、第1ストークスベクトルS、第2ストークスベクトルSのそれぞれに対して観測期間内における時間平均値<S>、<S>を算出し、それぞれの要素値である第0~第3ストークスパラメータの和を第0時間平均値<s10+s20>、第1時間平均値<s11+s21>、第2時間平均値<s12+s22>、第3時間平均値<s13+s23>として算出する(図18、ステップS124)。
 そして、変動特性解析部218は、(11)式に示すように第1時間平均値、第2時間平均値および第3時間平均値に対する二乗和の平方根を、第0時間平均値で除算して得られる値をTPUとして定める(図18、ステップS126)。TPUは、偏光状態の時間変動が少ないほど小さい値をとる。そのため、試料とする生体組織の活動が活発なほどTPUが小さくなることが期待される。これに対し、CPR、LPRまたは複屈折に対する上記の時間変動特性値は、偏光状態の時間変動が著しいほど大きい値をとる。
Figure JPOXMLDOC01-appb-M000011
 なお、測定ジョーンズ行列には、(12)式に示すようにノイズ成分が含まれる。(12)式において、n1H,n1V,n2H,n2Vは、それぞれ信号成分E1H,E1V,E2H,E2Vに加わるノイズ成分を示す。
Figure JPOXMLDOC01-appb-M000012
 偏光解析部216は、(13)式に示すように第1ストークスベクトルの時間平均値<S>、第2ストークスベクトルの時間平均値<S>について、それぞれに含まれる要素値の時間平均パワーからノイズ成分の要素値の時間平均パワーを差し引くことによりノイズ成分を補償することができる。ここで、変動特性解析部218は、(11)式において、ノイズ成分の補正前の第1ストークスベクトルの時間平均値<S>、第2ストークスベクトルの時間平均値<S>に代え、ノイズ成分の補正後の信号成分に対する第1ストークスベクトルの時間平均値<S’>、第2ストークスベクトルの時間平均値<S’>を代入してノイズ成分が除去されたTPUを算出することができる。
Figure JPOXMLDOC01-appb-M000013
 その他、変動特性解析部218は、所定の観測期間内における時間変動特性値として、例えば、観測点において定めたジョーンズ行列のフォンノイマンエントロピー(von Neumann Entropy)を算出してもよい。フォンノイマンエントロピーは、(14)式により定義される。(14)式において、λは、エルミート行列T(後述)の固有値を示し、個々の固有値λをそれらの固有値の総和で正規化して正規化固有値λ’が得られる。
Figure JPOXMLDOC01-appb-M000014
Figure JPOXMLDOC01-appb-M000015
 (15)式に示すようにLPRの期待値E(R)は個々の固有値に対応する光軸のLPRであるRの加重平均値となることが知られている。加重平均における個々の光軸のLPRに乗算される重み係数が、その光軸に対応する正規化固有値λ’で与えられる。そのため、正規化固有値λ’を用いて定義されるジョーンズ行列のフォンノイマンエントロピーも偏光状態に係る時間変動特性値の一種とみることができる。ジョーンズ行列のフォンノイマンエントロピーは、次に示す文献に詳しく記載されているが、本実施形態では、偏光特性値としてジョーンズ行列の時間変動による乱雑性(Randomness)を評価する。
 Masahiro Yamanari et al.: “Estimation of Jones matrix, birefringence and entropy using Cloude-Pottier decomposition in polarization-sensitive optical coherence tomography,” Biomedical Optics EXPRESS Vol. 7, No.9, p.3551-3572, published 1 Sep. 2016
 Masahiro Yamanari et al.: “Estimation of Jones matrix, birefringence and entropy using Cloude-Pottier decomposition in polarization-sensitive optical coherence tomography: erratum,” Biomedical Optics EXPRESS Vol. 7, No.11, p.4636-4637, published 1 Nov. 2016
 偏光解析部216と変動特性解析部218は、次の手順でフォンノイマンエントロピーHを算出することができる。まず、偏光解析部216は、観測点ごとに計測ジョーンズ行列Jmzの2行2列の要素を各行各列の順に並び替えて構成される4次元のベクトル[g1H g2H g1V g2Vを目標ベクトルκとして構成する。この目標ベクトルκは、計測ジョーンズ行列Jmzと実質的に同じ値を示す。
 偏光解析部216は、目標ベクトルのエルミート共役κ にその目標ベクトルκを乗じて4行4列の正方行列κκ を算出する。変動特性解析部218は、(16)式に示すように所定期間内における正方行列κκ の時間平均値<κκ >を行列Tとして定める。行列Tは、4行4列の半正定値エルミート行列(Positive Semidefinite Hermitian Matrix)となる。
Figure JPOXMLDOC01-appb-M000016
 そして、変動特性解析部218は、行列Tを対角化し、4個の固有値λ(i=1~4)を算出する。但し、iは、固有値λの降順に定められるインデックスである。変動特性解析部218は、個々の固有値λを4個の固有値の総和Σj=1 λで除算して正規化固有値λ’を算出する。変動特性解析部218は、(14)式に示すように、正規化固有値λ’とその対数値log(λ’)との積の総和を正負反転して得られる値をフォンノイマンエントロピーHとして算出することができる。フォンノイマンエントロピーは0以上1以下の値をとる。完全にランダムなジョーンズ行列のセットに対しては、フォンノイマンエントロピーHは1となる。但し、(14)式において対数値の底を4とする。
 なお、変動特性解析部218は、(17)式に示すように、測定ジョーンズ行列に基づいて定めたフォンノイマンエントロピーHから、ノイズ成分のフォンノイマンエントロピーHを差し引いて、信号成分のフォンノイマンエントロピーHを定めてもよい。
Figure JPOXMLDOC01-appb-M000017
 変動特性解析部218は、(18)式に示すように、ノイズ成分のフォンノイマンエントロピーH(E,E)をジョーンズ行列Jmzの部分空間をなす第1ジョーンズベクトルE、第2ジョーンズベクトルEのそれぞれに加わる第1ノイズ成分nのエントロピーである第1ノイズ成分エントロピーH(E)と第2ノイズ成分nのエントロピーである第2ノイズ成分エントロピーH(E)の和に近似することができる。但し、第1ノイズ成分nと第2ノイズ成分nは互いに独立と仮定する。
Figure JPOXMLDOC01-appb-M000018
 変動特性解析部218は、(19)式に示すように、第iノイズ成分エントロピーH(E)を、第iノイズ成分の第j固有値ζ (i)とその対数値log(ζ (i))との総和を正負反転して算出することができる。
Figure JPOXMLDOC01-appb-M000019
 但し、変動特性解析部218は、(20)式に示すように、1に第iジョーンズベクトルのTPUであるP(i)を加算もしくは減算した値に2で除算した値をそれぞれ第1固有値ζ (i)、第2固有値ζ (i)として算出することができる。
Figure JPOXMLDOC01-appb-M000020
 変動特性解析部218は、(21)式に示すように、第iジョーンズベクトルのTPUであるP(i)を、上記のノイズ成分で補正後の第1ストークスベクトル、第2ストークスベクトルのそれぞれの第1~第3ストークスパラメータの時間平均値の平方和の平方根を第0ストークスパラメータの時間平均値で除算して算出することができる。
Figure JPOXMLDOC01-appb-M000021
(時間変動特性値の算出例)
 次に、上記の時間変動特性値の算出例について説明する。図13は、複屈折率の分散の空間分布の一例を示す図である。図13では、生体組織内の観測点ごとの複屈折率の分散を濃淡で示す。観測点は、x-z平面内に分布している。明るい部分ほど複屈折率の分散が大きく、暗い部分ほど複屈折率の分散が小さい。図13において上方の黒で塗りつぶされている部分は組織外を示す。全体として組織の表面よりも内部の方が複屈折率の分散が大きくなる傾向がある。
 図14は、複屈折率の分散の空間分布の他の例を示す図である。図14は、生体組織内の観測点ごとの複屈折率の分散を濃淡で示す。図15は、観測点ごとの平均局所複屈折率(mean local birefringence)を示す。図15に示す例では、図14と同じ生体組織を観測対象としている。但し、図14に示す例では、図13とは異なる生体組織を観測対象としている。図14、図15ともに周囲よりも明るい部分が生体組織の分布範囲を示す。いずれも、図面の下方において複屈折率または平均局所複屈折率が相対的に高く、図面の左方において複屈折率または平均局所複屈折率が相対的に高い傾向がある。図16は、複屈折率の分散と平均局所複屈折率との相関性の一例を示す図である。図16の縦軸、横軸にそれぞれ複屈折率の分散、平均局所複屈折率を示す。図16は、局所複屈折率と複屈折率の分散との間に有意な相関性があることを示す。相関係数は、0.776となった。このことは、局所複屈折率が高い部分ほど複屈折率の分散が大きい傾向があることを裏付ける。図17は、複屈折率の分散と対数強度の分散との相関性の一例を示す図である。図17の縦軸、横軸にそれぞれ複屈折率の分散、対数強度の分散を示す。対数強度は、観測点ごとの信号強度の対数値である。図17は、局所複屈折率と複屈折率の分散との間に有意な相関性がないことを示す。相関係数は、0.280となった。このことは、対数強度だけでは偏光状態を表現しきれないことを裏付ける。
 図19は、TPUの空間分布の一例を示す図である。図19は、生体組織内の観測点ごとのTPUを濃淡で示す。図19に示す例では、図13と同じ生体組織を観測対象としている。明るい部分ほどTPUが大きく、暗い部分ほどTPUが小さいことを示す。全体として組織の表面よりも内部においてTPUが小さくなる傾向がある。この傾向は、複屈折率の分散とは逆の傾向を示す。TPUが大きいほど偏光状態の時間変化が少ないことを示すのに対し、複屈折率の分散が大きいほど偏光状態の時間変化が大きいことによる。
 図20は、複屈折率の分散の時間変化の一例を示す図である。図21は、TPUの時間変化の他の例を示す図である。図22は、複屈折率の動的コントラストの時間変化の一例を示す図である。図23は、対数強度分散の時間変化の一例を示す図である。図20、図21、図22、図23は、共通の生体組織内の、ある観測点における1時間ごとの複屈折率の分散、TPU、複屈折率、対数強度分散をそれぞれ示す。当該生体組織は時間経過により活動状態が低下する。
 図20は、時刻0から44時間経過後までは、時間経過により複屈折率の分散が有意に減少する傾向を示すのに対し、45時間経過後から60時間経過後までは複屈折率がほぼ一定となる。図20は、時刻0から44時間経過後までは、時間経過により複屈折率の分散が有意に減少する傾向を示す。相関係数は、-0.9486となった。これに対し、45時間経過後から60時間経過後までは複屈折率がほぼ一定となった。相関係数は、-0.1711となった。図21は、時刻0から42時間経過後までは、時間経過によりTPUが有意に増加する傾向を示す。相関係数は、0.9413となった。これに対し、43時間経過後から60時間経過後まではTPUがほぼ一定となった。相関係数は、0.0735となった。図22は、時間経過により複屈折率の動的コントラストが有意に減少する傾向を示す。相関係数は、-0.905となった。図20は、時刻0から3時間経過後までは対数強度分散が減少し、4時間経過後から19時間経過後までは対数強度分散が増加し、20時間経過後から42時間経過後までは対数強度分散が減少し、42時間経過後から60時間経過後までは対数強度分散がほぼ一定となった。これらの例は、本実施形態によれば試料内の観測点ごとに偏光特性値の時間変動特性値を計測することができ、計測された時間変動特性値に基づいて生体組織の活動状態を評価することができることを示す。
 以上に説明したように、上述した実施形態に係る計測信号処理装置200は、試料に第1偏光状態で入射される第1入射成分(例えば、水平偏光成分)に対して試料から反射または散乱した成分に第1入射成分が干渉した第1干渉成分のうち第1偏光状態(例えば、水平偏光)を有する第1計測信号(例えば、第1水平偏波スペクトル干渉信号)、第1干渉成分に対する第2偏光状態(例えば、垂直偏光)を有する第2計測信号(例えば、第2水平偏波スペクトル干渉信号)、試料に第2偏光状態で入射される第2入射成分(例えば、垂直偏光成分)に対して試料から反射または散乱した成分に第1入射成分が干渉した第1干渉成分のうち第1偏光状態を有する第3計測信号(例えば、第1垂直偏波スペクトル干渉信号)、および第2干渉成分に対する第2偏光状態を有する第4計測信号(例えば、第2垂直偏波スペクトル干渉信号)に基づいて、試料内の観測点における偏光特性に基づく偏光特性値を定める偏光解析部216と、偏光特性値の時間変動特性を示す時間変動特性値を定める変動特性解析部218と、を備える。この構成により、試料内の観測点ごとに偏光特性値の時間変動特性値が定まる。偏光特性値の時間変動は組織の活性と有意な相関があるため、観測点ごとに定まった時間変動特性値の分布により組織の活性を定量的に評価することができる。
 また、偏光解析部216は、第1計測信号、第2計測信号、第3計測信号および第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、試料内の観測点におけるジョーンズ行列と試料の表面におけるジョーンズ行列から観測点における累積ジョーンズ行列を定め、偏光特性値として、累積ジョーンズ行列の固有値間の位相差であるCPRを定めてもよい。この構成により、試料の表面と観測点の間の偏光成分間の位相差に相当するCPRの時間変動特性値が組織の活性の評価に利用することができる。
 また、偏光解析部216は、第1計測信号、第2計測信号、第3計測信号、および第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、試料内の第1観測点におけるジョーンズ行列と試料内の第2観測点におけるジョーンズ行列から第1観測点ならびに第2観測点との間の局所ジョーンズ行列を定め、局所ジョーンズ行列の2個の固有値間の位相差であるLPRを定めてもよい。この構成により、第1観測点と第2観測点の区間に生じる偏光成分間の位相差に相当するLPRの時間変動特性値が組織の活性の評価に利用することができる。そのため、微小な領域ごとに組織の活性を評価することができる。
 また、偏光解析部216は、局所位相遅延を試料に入射する入射光の波数と、第1観測点と第2観測点の厚みで除算して複屈折率を定めてもよい。そのため、微小な領域ごとに組織の活性を評価することができるとともに、複屈折率を用いることで観測対象とする他の組織または既存の物質との複屈折率の比較が容易となる。
 変動特性解析部218は、偏光特性値の分散または標準偏差に基づく時間変動特性値を算出してもよい。そのため、偏光特性値の時間変動の大きさが定量的に評価される。
 変動特性解析部218は、偏光特性値の対数値の分散または標準偏差に基づく時間変動特性値を算出してもよい。分散または標準偏差の算出の過程で、偏光特性値に潜在的に乗じられた定数が消去されるため、偏光特性値の実質的な時間変動が評価される。また、対数値をとることで、観測対象とする他の組織または既存の物質とのスケールの異なる時間変動特性値の比較が容易になる。
 変動特性解析部218は、偏光特性値の標準偏差を複屈折率の平均値で除算して動的コントラストを算出してもよい。偏光特性値の平均値で除算することで、偏光特性値の標準偏差が正規化されるので、スケールの変化を伴わずに偏光特性値の実質的な時間変動が評価される。
 偏光解析部216は、偏光特性値として、前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルを、それぞれ第1ストークスベクトル、および第2ストークスベクトルに変換し、変動特性解析部218は、時間変動特性値として、第1ストークスベクトルの時間平均値と前記第2ストークスベクトルの時間平均値に基づいてTPUを定めてもよい。この構成によれば、TPUを用いることで観測点における偏光状態の時間経過による一様性が定量化される。TPUは、組織の活性が低いほど大きくなる傾向がある。そのため、観測点ごとに定まったTPUの分布により組織の不活発な状態を定量的に評価することができる。
 変動特性解析部218は、第1ストークスベクトルからノイズ成分を差し引いた補正後の第1ストークスベクトルの時間平均値と、第2ストークスベクトルからノイズ成分を差し引いた補正後の第2ストークスベクトルの時間平均値に基づいてTPUを定めてもよい。この構成によれば、第1ストークスベクトルと第2ストークスベクトルからノイズ成分が除去され、信号成分が残される。そのため、TPUに対するノイズの影響を抑制して、組織の活性の評価を的確に行うことができる。
 偏光解析部216は、前記偏光特性値として、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号から観測点ごとにジョーンズ行列を定め、変動特性解析部218は、前記時間変動特性値として、前記ジョーンズ行列のフォンノイマンエントロピーを算出してもよい。この構成によれば、観測点における偏光状態を示すジョーンズ行列の乱雑性が定量化される。そのため、観測点ごとに定まったフォンノイマンエントロピーの分布により組織の活性を定量的に評価することができる。
 変動特性解析部218は、前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルから、それぞれ変換された第1ストークスベクトルの時間偏光均一度、および第2ストークスベクトルの時間偏光均一度からノイズ成分のエントロピーを算出し、フォンノイマンエントロピーをノイズ成分のエントロピーに基づいて補正してもよい。この構成によれば、フォンノイマンエントロピーからノイズ成分のエントロピーの寄与が補償され、信号成分のフォンノイマンエントロピーが得られる。そのため、フォンノイマンエントロピーに対するノイズの影響を抑制して、組織の活性の評価を的確に行うことができる。
 計測信号処理装置200は、入力値の変化に対して単調に変化する出力値を与える関数を用いて観測点ごとの時間変動特性値に対する出力値を信号値として有する画像データを生成する画像処理部220を備えてもよい。この構成によれば、観測点ごとの時間変動特性値に対応する輝度または色調の分布を有する画像が得られる。計測信号処理装置200は、時間変動特性値に基づいて試料とする生体組織の活性状態を示す評価値を定める出力処理部222を備えてもよい。かかる評価値は、例えば、活性度が高いほど大きい実数値であればよい。出力処理部222には、例えば、予め評価値と時間変動特性値と関係を示す関数およびそのパラメータを設定しておく。出力処理部222は、その評価値を記憶部230に記憶してもよいし、他の機器に出力してもよい。画像処理部220は、出力処理部222が観測点ごとに算出した評価値を上記のように画素ごとの信号値に変換し、変換した信号値を有する画像データを生成してもよい。そのため、利用者は得られた画像に接することで観測領域における組織の活性を容易に評価することができる。
 以上、図面を参照してこの発明の実施形態について詳しく説明してきたが、具体的な構成は上述のものに限られることはなく、この発明の要旨を逸脱しない範囲内において様々な設計変更等をすることが可能である。
 例えば、上記の説明では、評価装置20、計測信号処理装置200が、それぞれOCTシステム1の一部である場合を例にしたが、これには限られない。評価装置20、計測信号処理装置200は、それぞれOCTシステム1から独立し、光学系を備えない単一の機器であってもよい。その場合、計測信号処理装置200の制御部210において光学系制御部212が省略されてもよい。同様に、評価装置20には、光学系を制御するための制御手段が省略されてもよい。評価装置20、計測信号取得部214は、それぞれ光学系に限られず、データ蓄積装置、PCなど、他の機器から検出信号と計測信号を無線または有線で、例えば、ネットワークを経由して取得してもよい。
 また、計測信号処理装置200は、上記のように入出力部240、表示部250および操作部260のいずれも備えてもよいし、それらの一部または全部が省略されてもよい。また、評価装置20も、入出力部240、表示部250および操作部260に相当する機能構成のいずれも備えてもよいし、それらの一部または全部が省略されてもよい。
 また、計測信号処理装置200の制御部210において画像処理部220と出力処理部222の一方または両方が省略されてもよい。画像処理部220が省略される場合、制御部210は、生成した時間変動特性値を示すデータを、データ蓄積装置、PC、他の画像処理装置など、他の機器に無線または有線で、例えば、ネットワークを経由して出力してもよい。また、評価装置20も、画像処理部220と出力処理部222の一方または両方に相当する機能構成のいずれも備えてもよいし、それらの一部または全部が省略されてもよい。出力先とする機器は、画像処理部220と同様の機能、つまり、評価装置20または計測信号処理装置200から入力されるデータに基づいて出力画像データを生成し、生成した出力画像データに基づく画像を表示する機能を有してもよい。評価装置20において、計測信号処理装置200の画像処理部220と出力処理部222の一方または両方に相当する機能構成が省略されてもよい。
 また、上述した実施形態における評価装置20もしくは計測信号処理装置200の一部、または全部を、LSI(Large Scale Integration)等の集積回路として実現してもよい。評価装置20もしくは計測信号処理装置200の各機能ブロックは個別にプロセッサ化してもよいし、一部、または全部を集積してプロセッサ化してもよい。また、集積回路化の手法はLSIに限らず専用回路、または汎用プロセッサで実現してもよい。また、半導体技術の進歩によりLSIに代替する集積回路化の技術が出現した場合、当該技術による集積回路を用いてもよい。
 以上、本発明の好ましい実施形態を説明したが、本発明はこれら実施形態及びその変形例に限定されることはない。本発明の主旨を逸脱しない範囲で、構成の付加、省略、置換、およびその他の変更が可能である。
 また、本発明は前述した説明によって限定されることはなく、添付の特許請求の範囲によってのみ限定される。
 上記の実施形態によれば、例えば、再生医療用培養組織、オルガノイドの品質管理、動物実験、培養組織による薬効測定の際の薬効評価、眼底視細胞の遺伝子治療における治療効果評価等において、きわめて有用である。
1…OCTシステム、10…撮影手段、20…評価装置、22…計測手段、24…評価手段、102…光源、110…偏光遅延ユニット、128…プローブ、130…参照アーム、150…偏光分離検出ユニット、200…計測信号処理装置、210…制御部、212…光学系制御部、214…計測信号取得部、216…偏光解析部、218…変動特性解析部、220…画像処理部、222…出力処理部、230…記憶部、240…入出力部、250…表示部、260…操作部

Claims (23)

  1.  試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、
     前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測部と、
     前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価部と、
     を有する評価装置。
  2.  前記評価部は、
     前記時間変動特性値として、前記信号値の分散を算出する
     請求項1に記載の評価装置。
  3.  前記評価部は、前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差の二乗和を前記所定期間におけるフレーム数で除算して前記分散を観測点ごとに算出する
     請求項2に記載の評価装置。
  4.  前記評価部は、
     前記信号値と、当該信号値を時間シフト量τで時間シフトした信号値との相関係数を前記時間シフト量τごとに算出し、前記時間シフト量τの増加に伴う前記相関係数の減衰速度を前記時間変動特性値として算出する
     請求項1に記載の評価装置。
  5.  前記評価部は、
     前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差の二乗和を分散として算出し、
     前記所定期間内の各フレーム時刻における前記OCT信号の信号強度と当該信号強度の平均値との偏差と、当該フレーム時刻から時間シフト量τで時間シフトしたシフト時刻における前記OCT信号の信号強度と、時間シフトした当該信号強度の平均値との偏差との積の総和を共分散として算出し、
     前記共分散を前記分散で除算して前記相関係数として、前記時間シフト量τごとに算出し、
     前記時間シフト量τごとの前記相関係数を用いて所定の減衰関数を用いて回帰分析し、前記相関係数を近似する当該減衰関数のパラメータを前記減衰速度として観測点ごとに算出する
     請求項4に記載の評価装置。
  6.  前記評価部は、
     前記時間シフト量τを非零として算出される前記相関係数を用いて、前記減衰速度を算出する
     請求項4または請求項5に記載の評価装置。
  7.  前記計測部は、
     前記試料に第1偏光状態で入射される第1入射成分に対して前記試料から反射または散乱した成分に前記第1入射成分が干渉した第1干渉成分のうち第1偏光状態を有する第1計測信号、前記第1干渉成分に対する第2偏光状態を有する第2計測信号、前記試料に第2偏光状態で入射される第2入射成分に対して前記試料から反射または散乱した成分に前記第2入射成分が干渉した第2干渉成分のうち第1偏光状態を有する第3計測信号、および前記第2干渉成分に対する第2偏光状態を有する第4計測信号に基づいて、前記試料内の観測点における偏光特性に基づく偏光特性値を定め、
     前記評価部は、
     前記偏光特性値の時間変動特性を示す前記時間変動特性値を定める
     請求項1に記載の評価装置。
  8.  前記計測部は、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、前記試料内の観測点におけるジョーンズ行列と前記試料の表面におけるジョーンズ行列から前記観測点における累積ジョーンズ行列を定め、
     前記偏光特性値として、前記累積ジョーンズ行列の固有値間の位相差である累積位相遅指標値を定める
     請求項7に記載の評価装置。
  9.  前記計測部は、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、前記試料内の第1観測点におけるジョーンズ行列と前記試料内の第2観測点におけるジョーンズ行列から前記第1観測点ならびに前記第2観測点との間の局所ジョーンズ行列を定め、
     前記局所ジョーンズ行列の固有値間の位相差である局所位相遅延に基づく前記偏光特性値を定める
     請求項7に記載の評価装置。
  10.  前記計測部は、
     前記局所位相遅延を前記試料に入射する入射光の波数と、前記第1観測点と前記第2観測点との厚みで除算して複屈折率を定める
     請求項9に記載の評価装置。
  11.  前記評価部は、前記偏光特性値の分散または標準偏差に基づく前記時間変動特性値を算出する
     請求項10に記載の評価装置。
  12.  前記評価部は、前記偏光特性値の対数値の分散または標準偏差に基づく前記時間変動特性値を算出する
     請求項11に記載の評価装置。
  13.  前記評価部は、前記偏光特性値の標準偏差を前記複屈折率の平均値で除算して動的コントラストを算出する
     請求項11に記載の評価装置。
  14.  前記計測部は、前記偏光特性値として、前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルを、それぞれ第1ストークスベクトル、および第2ストークスベクトルに変換し、
     前記評価部は、前記時間変動特性値として、前記第1ストークスベクトルの時間平均値と前記第2ストークスベクトルの時間平均値に基づいて時間偏光均一度を定める
     請求項7に記載の評価装置。
  15.  前記計測部は、前記第1ストークスベクトルからノイズ成分を差し引いた補正後の第1ストークスベクトルの時間平均値と、前記第2ストークスベクトルからノイズ成分を差し引いた補正後の第2ストークスベクトルの時間平均値に基づいて時間偏光均一度を定める
     請求項14に記載の評価装置。
  16.  前記計測部は、前記偏光特性値として、前記第1計測信号、前記第2計測信号、前記第3計測信号、および前記第4計測信号に基づいて、観測点ごとにジョーンズ行列を定め、
     前記評価部は、前記時間変動特性値として、前記ジョーンズ行列のフォンノイマンエントロピーを算出する
     請求項7に記載の評価装置。
  17.  前記評価部は、
     前記第1計測信号と前記第2計測信号に基づく第1ジョーンズベクトル、および前記第3計測信号と前記第4計測信号に基づく第2ジョーンズベクトルから、それぞれ変換された第1ストークスベクトルの時間偏光均一度、および第2ストークスベクトルの時間偏光均一度からノイズ成分のエントロピーを算出し、
     前記フォンノイマンエントロピーを前記ノイズ成分のエントロピーに基づいて補正する
     請求項16に記載の評価装置。
  18.  前記第1偏光状態は水平偏光であり、前記第2偏光状態は垂直偏光であり、
     前記第1計測信号は第1水平偏波スペクトル干渉信号であり、
     前記第2計測信号は第2水平偏波スペクトル干渉信号であり、
     前記第3計測信号は第1垂直偏波スペクトル干渉信号であり、
     前記第4計測信号は第2垂直偏波スペクトル干渉信号である
     請求項7から請求項17のいずれか一項に記載の評価装置。
  19.  前記評価部は、
     前記所定期間よりも長い観測期間間隔ごとに前記時間変動特性値を算出する
     請求項1から請求項18のいずれか一項に記載の評価装置。
  20.  前記時間変動特性値に基づいて前記試料の活性状態を示す評価値を定める出力処理部を備える
     請求項1から請求項19のいずれか一項に記載の評価装置。
  21.  入力値の変化に対して単調に変化する出力値を与える関数を用いて前記観測点ごとの前記時間変動特性値に対する出力値を信号値として有する画像データを生成する画像処理部を備える
     請求項1から請求項20のいずれか一項に記載の評価装置。
  22.  評価装置における方法であって、
     試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、
     前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測ステップと、
     前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価ステップと、
     を有する評価方法。
  23.  評価装置のコンピュータに、
     試料とする生体組織の状態を示す光干渉断層(OCT:Optical Coherence Tomography)信号を取得し、
     前記OCT信号に基づく信号値を、前記試料における観測点ごとに取得する計測手順と、
     前記信号値の所定期間内における時間変動特性を示す時間変動特性値を算出する評価手順と、
     を実行させるためのプログラム。
PCT/JP2020/042497 2019-11-15 2020-11-13 評価装置、評価方法およびプログラム WO2021095868A1 (ja)

Priority Applications (2)

Application Number Priority Date Filing Date Title
US17/776,191 US20220390357A1 (en) 2019-11-15 2020-11-13 Evaluation device, evaluation method, and program
JP2021556189A JPWO2021095868A1 (ja) 2019-11-15 2020-11-13

Applications Claiming Priority (4)

Application Number Priority Date Filing Date Title
JP2019207348 2019-11-15
JP2019-207348 2019-11-15
JP2020-070309 2020-04-09
JP2020070309 2020-04-09

Publications (1)

Publication Number Publication Date
WO2021095868A1 true WO2021095868A1 (ja) 2021-05-20

Family

ID=75912772

Family Applications (1)

Application Number Title Priority Date Filing Date
PCT/JP2020/042497 WO2021095868A1 (ja) 2019-11-15 2020-11-13 評価装置、評価方法およびプログラム

Country Status (3)

Country Link
US (1) US20220390357A1 (ja)
JP (1) JPWO2021095868A1 (ja)
WO (1) WO2021095868A1 (ja)

Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2019080724A (ja) * 2017-10-30 2019-05-30 キヤノン株式会社 画像処理装置、画像処理方法及びプログラム
JP2019170706A (ja) * 2018-03-28 2019-10-10 株式会社トプコン 眼科装置、及び眼科情報処理プログラム
JP2019170997A (ja) * 2018-03-29 2019-10-10 イメドース システムズ ゲーエムベーハー 網膜の血管の内皮機能を検査するための方法及び装置

Patent Citations (3)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP2019080724A (ja) * 2017-10-30 2019-05-30 キヤノン株式会社 画像処理装置、画像処理方法及びプログラム
JP2019170706A (ja) * 2018-03-28 2019-10-10 株式会社トプコン 眼科装置、及び眼科情報処理プログラム
JP2019170997A (ja) * 2018-03-29 2019-10-10 イメドース システムズ ゲーエムベーハー 網膜の血管の内皮機能を検査するための方法及び装置

Non-Patent Citations (2)

* Cited by examiner, † Cited by third party
Title
APELIAN, C. ET AL.: "Dynamic full field optical coherence tomography: subcellular metabolic contrast revealed in tissues by interferometric signals temporal analysis", BIOMEDICAL OPTICS EXPRESS, vol. 7, no. 4, 24 March 2016 (2016-03-24), pages 1511 - 1524, XP055805619, DOI: https://doi.org/10.1364/BOE.7.001511 *
HUANG, X. R. ET AL.: "Temporal change of retinal nerve fiber layer reflectance speckle in normal and hypertensive retinas", EXPERIMENTAL EYE RESEARCH, vol. 186, no. 107738, 17 July 2019 (2019-07-17), pages 1 - 11, XP085776185, DOI: https://doi.org/10.1016/j.exer.2019.107738 *

Also Published As

Publication number Publication date
JPWO2021095868A1 (ja) 2021-05-20
US20220390357A1 (en) 2022-12-08

Similar Documents

Publication Publication Date Title
US8208996B2 (en) Imaging of polarization scrambling tissue
US9471277B2 (en) Program for correcting data measured by PS-OCT and PS-OCT system equipped with the program
US10194803B2 (en) Control apparatus, measurement apparatus, control method, and storage medium
US7612880B2 (en) Advanced polarization imaging method, apparatus, and computer program product for retinal imaging, liquid crystal testing, active remote sensing, and other applications
KR101332222B1 (ko) 광간섭 단층촬영법을 이용해서 샘플 내에서 적어도 하나의 위치를 결정하는 방법, 시스템 및 그 방법을 구현하기 위한 소프트웨어가 저장되어 컴퓨터로 판독 가능한 매체
USRE46412E1 (en) Methods and systems for performing angle-resolved Fourier-domain optical coherence tomography
WO2020155415A1 (zh) 基于特征空间的光学相干层析的三维血流造影方法及系统
US8280494B2 (en) Apparatus and method to measure a spectroscopic characteristic in an object
US9226655B2 (en) Image processing apparatus and image processing method
JP6256879B2 (ja) 偏光感受型光画像計測システム及び該システムに搭載されたプログラム
US20180003479A1 (en) Image processing apparatus and image processing method
JP2008519264A (ja) 偏光感応性光コヒーレンストモグラフィを用いて偏光非解消の偏光パラメータを測定するジョーンズ行列に基づく解析を行うシステム及び方法
JP2010539491A (ja) 低コヒーレンス干渉法(lci)のための装置、システムおよび方法
US20140293289A1 (en) Method for Generating Two-Dimensional Images From Three-Dimensional Optical Coherence Tomography Interferogram Data
US9918623B2 (en) Optical tomographic imaging apparatus
KR101053222B1 (ko) 멀티라인 카메라를 이용한 광간섭성 단층촬영장치
US20170231484A1 (en) Image processing apparatus, image processing method, and storage medium
JP6648891B2 (ja) 物質含有量を断層可視化する装置および方法
JP2010151684A (ja) 局所的な複屈折情報を抽出可能な偏光感受光画像計測装置
WO2021095868A1 (ja) 評価装置、評価方法およびプログラム
JP5905711B2 (ja) 光画像計測装置
Verrier et al. Influence of interfaces reflectivity for central thickness measurement of a contact lens by low coherence interferometry
Thompson Interpretation and medical application of laser biospeckle
CN116849601A (zh) 一种拓展系统动态测量范围的方法、装置及电子设备
EP4447785A1 (en) Interferometer-based synthetic multi-exposure speckle imaging (symesi) method and system

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: 20887399

Country of ref document: EP

Kind code of ref document: A1

ENP Entry into the national phase

Ref document number: 2021556189

Country of ref document: JP

Kind code of ref document: A

NENP Non-entry into the national phase

Ref country code: DE

122 Ep: pct application non-entry in european phase

Ref document number: 20887399

Country of ref document: EP

Kind code of ref document: A1