Login| Sign Up| Help| Contact|

Patent Searching and Data


Title:
AUTOMATIC ASSESSMENT OF TIME-RESOLVED OCT IMAGES FOR SELECTIVE RETINA THERAPY
Document Type and Number:
WIPO Patent Application WO/2017/178059
Kind Code:
A1
Abstract:
The invention relates to a method for predicting a pattern in an optical coherence tomography sample image (OCT) data set (X), corresponding to a structural change in a retina (51) of an eye (5), comprising the steps of performing an OCT measurement of the retina (51), automatically providing a sample image data set (X) comprising a plurality of image data points (xi), automatically extracting a first plurality of features (u) from the sample image data set (X), and automatically predicting the presence or absence of the pattern from the first plurality of features (u) by means of a classification algorithm, wherein a second plurality of features (u) extracted from a training image data set serves as a training data set of the classification algorithm. The invention further relates to a system, and a computer program for predicting a pattern in an optical coherence tomography sample image data set (X).

Inventors:
SZNITMAN RAPHAEL (CH)
STEINER PATRICK (CH)
WOLF SEBASTIAN (CH)
ZBINDEN SARAH (CH)
Application Number:
PCT/EP2016/058295
Publication Date:
October 19, 2017
Filing Date:
April 14, 2016
Export Citation:
Click for automatic bibliography generation   Help
Assignee:
UNIVERSITÄT BERN (CH)
International Classes:
G06T7/00; A61B3/12
Other References:
ZBINDEN SARAH ET AL: "Automatic assessment of time-resolved OCT images for selective retina therapy", INTERNATIONAL JOURNAL OF COMPUTER ASSISTED RADIOLOGY AND SURGERY, SPRINGER, DE, vol. 11, no. 6, 11 April 2016 (2016-04-11), pages 863 - 871, XP035942132, ISSN: 1861-6410, [retrieved on 20160411], DOI: 10.1007/S11548-016-1383-6
STEINER ET AL.: "Time-Resolved Ultra-High Resolution Optical Coherence Tomography for Real-Time Monitoring of Selective Retina Therapy", INVESTIGATIVE OPHTHALMOLOGY & VISUAL SCIENCE, vol. 56, no. 11, pages 6654 - 6662
EBNETER, BERGER; ZINKERNAGEL; POVAZAY; MEIER; KOWAL; BRINKMANN; WOLF; SZNITMAN: "Time-Resolved Ultra-High Resolution Optical Coherence Tomography for Real-Time Monitoring of Selective Retina Therapy", INVESTIGATIVE OPHTHALMOLOGY AND VISUAL SCIENCE, vol. 56, no. 11, 2015, pages 6654 - 6662
STEINER ET AL.: "Retinal laser lesion visibility in simultaneous ultra-high axial resolution optical coherence tomography", IEEE PHOTONICS JOURNAL, vol. 6, no. 6, 2014, pages 1 - 11
Attorney, Agent or Firm:
SCHULZ, Ben Jesko (DE)
Download PDF:
Claims:
Claims

A method for predicting a pattern in an optical coherence tomography sample image data set (X), wherein the pattern corresponds to a structural change in a retina (51 ) of an eye (5), and wherein the method comprises the steps of

• performing an optical coherence tomography (OCT) measurement of said retina (51 ),

• automatically providing from said measurement a sample image data set (X) comprising a plurality of image data points (x,), wherein each image data point (x,) comprises a time coordinate value, a spatial coordinate value, and an intensity value,

• automatically extracting a first plurality of features (u) from said sample image data set (X),

• automatically predicting the presence or absence of said pattern in said sample image data set (X) from said first plurality of features (u) by means of a classification algorithm, wherein a second plurality of features (u) extracted from at least one training image data set serves as a training data set of said classification algorithm.

The method according to claim 1 , wherein a statistical measure,

particularly a mean, a median, a root mean square, a variance, or a standard deviation, or a function of said statistical measure, of the intensity values of at least a subset of said plurality of image data points (x,) is extracted as a feature (u) of said first plurality of features (u).

The method according to claim 1 or 2, wherein the intensity values of at least one subset of said plurality of image data points (x,) are aggregated into an aggregate value (a) by means of an aggregation function, particularly a root mean square function or a sum function, and wherein a statistical measure or a function of said statistical measure of said aggregate value (a) is extracted as a feature (u) of said first plurality of features (u).

The method according to claim 3, wherein the image data points (x,) of said subset of said plurality of image data points (x,) each comprise the same time coordinate value.

The method according to claim 3 or 4, wherein the number of aggregate values (a) exceeding a threshold subset value (sth) is extracted as a feature (u) of said first plurality of features (u).

6. The method according to one of the claims 3 to 5, wherein said aggregation function is described by the formula

St =∑Ll=1X{l, t),

wherein St is a group sum of intensity values X(l, t) calculated from groups of image data points (x,) having the same time coordinate value (t), wherein I labels the spatial coordinate value, and L is the maximum number of spatial coordinate values of said sample image data set (X), wherein particularly a respective offset term is subtracted from the respective group sum St, more particularly by a moving average operation, resulting in a plurality of corrected group sums (St,c), wherein the number of group sums (St) or the number of corrected group sums (St,c) exceeding a threshold value (sth) is extracted as a feature (u) of said first plurality of features.

7. The method according to one of the preceding claims, wherein the sample image data set (X) is split into a plurality of blocks (Bn) of image data points (x,), wherein the time coordinate values of the image data points (x,) in each block (Bn) are in a respective partial range, and wherein a respective block feature (un) of said first plurality of features (u) is extracted from at least one block (Bn).

8. The method according to claim 7, wherein the intensity values (x,) of at least one block (Bn) are aggregated into a respective block aggregate value (an) by means of an aggregation function, and wherein a statistical measure, particularly a mean, a root mean square, a variance, or a standard deviation, of the block aggregate value (an) is extracted as a block feature (un) of said first plurality of features (u).

9. The method according to claim 7 or 8, wherein a plurality of image data

points (xi 0) from said sample image data set (X) is allocated to a reference block (B0), and wherein for at least one image data point (x,) of each block (Bn) a corrected intensity value is generated by subtracting from the respective intensity value of the at least one image data point (x,) of each block (Bn) a reference intensity value of a corresponding image data point (xi 0) of the reference block (B0).

10. The method according to one of the claims 7 to 9, wherein at least one of said block features (un) is determined by means of the formula wherein I labels the spatial coordinate value, t labels the time coordinate value, L is the maximum number of spatial coordinate values of the respective block (Bn), Tn is the maximum number of time coordinate values of the respective block (Bn), Bn (l,t) is a reference-subtracted intensity value derived by subtraction of a respective intensity value taken from a reference block (B0) from the intensity value of the image data point (xi n) in the block Bn, and μη' is a mean determined by means of the formula

/½.' = r∑i tTn Bn' (I, t), and/ or wherein a maximum difference (ud) between subsequent block features (un, un+i) of subsequent blocks (Bn, Bn+1) is extracted as an additional feature (u) of said first plurality of features (u).

1 1. The method according to one of the claims 7 to 10, wherein at least one of said block features (un) is determined by means of the formula

wherein vn is a variance of the block Bn determined according to the formula

wherein v0 is a variance of a reference block (B0), and wherein I labels the spatial coordinate value, t and t' label the time coordinate values, Ln is the maximum number of spatial coordinate values of the respective block (Bn), Tn is the maximum number of time coordinate values of the respective block (Bn), Rn (t) is a root mean square of the intensity values Bn(l,t) of the image data points Rn(l,t) of the block Bn at the time coordinate value t, wherein Rn is determined by means of the formula

and/ or wherein a maximum value (umax) of the block features (un) is extracted as an additional feature of said first plurality of features (u).

12. The method according to one of the preceding claims, wherein at least one frequency value (f,), and/ or at least one phase value (p,) of at least one image data point (x,) is provided by means of

• a Fourier transformation, particularly a short time Fourier transformation, of said sample image data set (X), or • by means of a spectrometer (24), particularly comprised in a system (1 ) for predicting a pattern in a sample image data set (X) according to claim 14,

wherein at least one feature (u) of said first plurality of features (u) is extracted from said at least one frequency value (f,) and/ or said at least one phase value (pi).

13. The method according to claim 12, wherein

• a spectrogram (Z) comprising a plurality of time and frequency

dependent amplitude values is generated by means of a Fourier Transformation, particularly a Short-Time-Fourier-Transformation (ST FT) according to the formula Ζ(τ, ω) =∑=1 R(x)h(t - τ)θ_ίωτ, wherein ω is a frequency, t labels a first time coordinate value of the image data points (x,), T is the maximum number of time coordinate values of said sample image data set (X), h(t - τ) is a windowing function, T labels a second time coordinate value of the windowing function, R(x) is an aggregate value generated by means of a root mean square aggregation function, and wherein

• said spectrogram (Z) is converted to a sum spectrogram ( ) by means of the formula Ζ(τ) =∑ω Ζ(τ, ω), wherein particularly the amplitude values of said spectrogram (Z) are median filtered prior to conversion to said sum spectrogram ( ), and wherein

• said sum spectrogram is normalized, particularly by means of the

formula (τ) = zχW(Ζ(~τm)-inπά(zηW x (Ζ(τ))

' πια (Ζ)(τ)γ , wherein ma

J J is the maximum of said sum spectrogram (Z), and min(Z(x)) is the minimum of said sum spectrogram (Z), resulting in a normalized spectrogram (Z), and wherein

• a crest factor calculated according to the formula u = js

extracted as a feature (u) of said first plurality of features (u), wherein particularly the amplitude values of the spectrogram are converted to a dB scale prior to generating the crest factor,

and/ or wherein the sum (τ) is extracted as an additional feature (u) of said first plurality of features (u).

14. A system (1 ) for predicting a pattern in an optical coherence tomography

sample image data set (X), comprising i. an optical coherence tomography (OCT) device (2) comprising a first laser source (21 ), a measuring arm (22), a reference arm (23), which are optically connected, particularly by means of a fiber coupler (25), ii. an ophthalmoscope (27), particularly a scanning digital ophthalmoscope (SDO), which is coupled into the light path of said measuring arm (22) of said optical coherence tomography (OCT) device (2), such that a plurality of light intensities of reflected light, and/or diffracted light, and/or a plurality of phases in a retina (51 ) of an eye (5) observed by said ophthalmoscope (27) are measurable by means of said optical coherence tomography device (2),

iii. a second laser source (3), which is coupled into the light path of said measuring arm (22) of said optical coherence tomography device (2), particularly by means of a dichroic mirror (28), such that a retina (51 ) of an eye (5) is observable by means of said optical coherence tomography device (2) during illumination of said retina (51 ) by said second laser source (3)

characterized by

an analysis device (4), which is adapted to

• automatically extract a first plurality of features (u) from said sample image data set (X) received from said optical coherence tomography (OCT) device (2),

• automatically predict the presence or absence of a pattern in said

sample image data set (X) from said first plurality of features (u) by means of a classification algorithm, wherein said pattern corresponds to a structural change in the retina (51 ) of the eye (5), and wherein a second plurality of features (u) extracted from at least one training image data set serves as a training data set of said classification algorithm,

particularly by means of steps of the method according to one of the claims 1 to 13.

15. A computer program comprising program code, which is adapted to perform the steps of • automatically extracting a first plurality of features (u) from a sample image data set (X) received from an optical coherence tomography (OCT) device (2),

• automatically predicting the presence or absence of a pattern in said sample image data set (X) from said first plurality of features (u) by means of a classification algorithm, wherein said pattern corresponds to a structural change in a retina (51 ) of an eye (5), and wherein a second plurality of features (u) extracted from at least one training image data set serves as a training data set of said classification algorithm, when the program is executed on a computer, particularly by means of steps of the method according to one of the claims 1 to 13.

Description:
Automatic assessment of time-resolved OCT images for selective retina therapy

The invention relates to a method, system, and computer program for predicting a pattern in an optical coherence tomography (OCT) data set, which corresponds to a structural change in a retina of an eye, particularly a change as a result of selective retina therapy (SRT). In particular, the method according to the invention can be used to analyze patient OCT data during SRT.

Retinal photocoagulation is used in the treatment of various retinal pathologies, such as retinal tears, ischemic retina, and macular diseases such as chorioretinopathia centralis serosa and diabetic edema. In this method, light is applied to a patient's retina in order to destroy retinal pigment epithelial cells, which reduces the amount of hypoxia and the VEGF load. However, conventional retinal photocoagulation has adverse effects. For example, side effects in treatment of macular diseases include laser scotoma and denaturation of the photoreceptor layer.

In selective retina therapy (SRT) the retina is subjected to laser pulses of defined pulse energies, which allows to selectively target retinal epithelial cells, leaving surrounding retinal layers intact. Monitoring the effects of SRT in the patient is desirable, particularly because necessary pulse energies are highly variable in individual patients due to different melanin concentrations. However, the effects of SRT are non-detectable by standard non-invasive methods such as fundus photography.

Recently, a novel method for monitoring the effects of SRT by optical coherence tomography (OCT) during treatment was described (Steiner et al., Time-Resolved Ultra-High Resolution Optical Coherence Tomography for Real-Time Monitoring of Selective Retina Therapy, Investigative Ophthalmology & Visual Science 56(1 1 ), 6654- 6662). Therein, an OCT device is used to provide axial scans (M-scans) of a patient retina at specific spots, which are irradiated by an SRT treatment laser. A visible change in the OCT signal was observed following the application of SRT laser pulses. The occurrence of the signal change was correlated with the presence of lesions in the retina as judged by fundus fluorescein angiography. However, the OCT data had to be analyzed by a human observer.

Therefore, the problem to be solved by the present invention is to provide an automated, robust and efficient method, system, and computer program for predicting a pattern in an OCT data set, which corresponds to a structural change in a retina of an eye, particularly a change as a result of selective retina therapy (SRT). This problem is solved by the subject matter of the independent method claim 1 , system claim 14, and computer program claim 15. Specific embodiments of the method are claimed by the dependent claims 2 to 13.

According to a first aspect of the invention, a method for predicting a pattern in an optical coherence tomography image data set is provided. The pattern in the image data set corresponds to a structural change in a retina of an eye that can be detected by means of an optical coherence tomography measurement of the retina.

In particular, the structural change in the retina may be a result of retina therapy, more particularly selective retina therapy, wherein structural changes in the retinal pigment epithelium (RPE) are induced by means of a treatment laser beam, which is targeted at the retina. Without wishing to be bound by theory, the structural changes in the retina following retinal therapy may comprise the formation of microbubbles, particularly by evaporation of water, more particularly in the RPE.

The method comprises a first step of performing an optical coherence

tomography (OCT) measurement of the retina.

In particular, this measurement comprises a scan in an axial direction (A-scan), more particularly a time series of axial scans (M-scan) of the retina. In particular, such a scan may be achieved by means of a digital scanning ophthalmoscope according to techniques described in the prior art (Ebneter, Berger, Zinkernagel, Povazay, Meier, Kowal, Brinkmann, Wolf, and Sznitman: Time-Resolved Ultra-High Resolution Optical Coherence Tomography for Real-Time Monitoring of Selective Retina Therapy;

Investigative Ophthalmology and Visual Science, 56, 11, 6654 - 6662, 2015).

Therein, a measurement laser beam is targeted at an area of the retina, and the focus of the measurement laser beam is changed in the axial direction by means of the digital scanning ophthalmoscope, such that the measurement laser beam is focused at different layers of the retina. The digital scanning ophthalmoscope is positioned in the light path of a measuring arm of an OCT device. Such OCT devices are

interferometers comprising a light source, particularly a laser light source, a measuring arm, a reference arm, and a detector, wherein the light emanating from the light source is guided to the measuring arm, is at least partially reflected and/ or diffracted by the sample, and constructively or destructively interferes with the light in the reference arm. The detector of such a device is adapted to measure the light intensity of the light beam resulting from interference of the light beams traveling through the measuring arm and the reference arm. In addition, the detector may comprise a spectrometer, which is adapted to split the light beam into light of different frequencies. In particular, the frequency distribution of the light beam reflected and/ or diffracted from the sample may be measured by means of the spectrometer. In particular, the phase shift between the light emanating from the light source and the light reflected and/ or diffracted by the sample may be measured by means of the spectrometer. Here, the sample refers to a portion of the retina, onto which the measuring light beam is focused by means of the ophthalmoscope.

This portion of the retina may reflect, diffract, and/ or absorb at least a part of the incoming light, depending on its optical properties, wherein the optical properties influence the amount, angle, phase shift, and/ or frequency distribution of the reflected, diffracted, and/ or absorbed light. In particular, optical properties include the refractive index and the absorption coefficient at a given frequency of light. In particular, these optical properties depend on the physical, chemical, and biological structure of the respective portion of the retina, wherein particularly the presence of microbubbles in the RPE may influence the optical properties.

In addition to the described optical properties, the OCT signal of a given portion of the retina measured at the detector at a given time also depends on the ratio of the length of the measuring arm and the length of the reference arm. In particular, the length of the measuring arm depends on the distance of the respective portion of the retina to the ophthalmoscope. In particularly, the length of the reference arm may be

changeable, particularly by means of an adjustable mirror.

In an M-scan, a time series of A-scans of the retina is performed, typically at a frequency of around 70 kHz, wherein each A-scan comprises measured intensity values of different portions of the retina corresponding to different axially arranged layers of the retina in the area targeted by the measuring laser beam. In particular, during an M-scan, the distance of the respective portion of the retina from the ophthalmoscope, and/ or the biological structure of the respective portion of the retina may change, particularly due to thermal fluctuations, or due to the effects of a retinal therapy, resulting in a temporal alteration of the OCT signal.

A second step of the method comprises automatically providing from the measurement a sample image data set comprising a plurality of image data points, wherein each image data point comprises a time coordinate value, a spatial coordinate value, and an intensity value.

Therein, each image data point corresponds to a respective portion of the retina at a given time, wherein the intensity value of the respective image data point is a result of the optical properties of the respective portion of the retina, and the ratio of the length of the measuring arm and the reference arm. The spatial coordinate corresponds to the axial position of the respective portion of the retina during an A-scan. The temporal coordinate corresponds to the time at which the respective A-scan was performed.

In particular, the image data points may further comprise at least one frequency value, and/ or at least one phase value. More particularly, a frequency value may comprise a signal amplitude at a given frequency of light, and a phase value may comprise a phase shift between the measuring laser beam and the reflected light.

In particular, the sample image data set is provided by transmitting data from the detector of the OCT device to an analysis device, particularly according to methods of the prior art. More particularly, the measured values may be converted to digital values by means of an analog digital converter.

In particular, the intensity values, time coordinate values, and/ or spatial coordinate values may be stored and/ or further processed by a data processing unit. More particularly, the intensity values may be stored and/ or displayed as grey values, for example 8 bit or 16 bit grey values. In particular, the image data set may be displayed as a plurality of pixels, wherein each pixel corresponds to a respective image data point, and wherein the grey value of each pixel corresponds to the intensity value of the respective image data point, and wherein the x-coordinate of the displayed image data set corresponds to the temporal coordinate value, and wherein the y-coordinate of the displayed image data set corresponds to the spatial coordinate value. Such a plurality of pixels is also termed "speckle pattern".

In the context of the present specification, the term "pattern" designates a specific property of the distribution of intensity values and/ or frequency values, and/ or phase values along the time coordinate and the spatial coordinate, and the term "speckle pattern" designates a specific distribution of intensity values along the time coordinate and the spatial coordinate. In particular, a "pattern" may be periodic or non-periodic with respect to the time coordinate.

The method comprises a third step of automatically extracting a first plurality of features from the sample image data set.

In the context of the present specification, the term "feature" refers to a mathematical value, particularly a real number or an integer, wherein the value reflects at least one property of an image data set. In particular, a "feature" is calculated by means of a mathematical function or a series of functions, wherein intensity values, frequency values, phase values, time coordinate values, and/ or spatial coordinate values of at least one image data point of the image data set are used as an input of the function or the series of functions.

The method comprises a fourth step of automatically predicting the presence or absence of the pattern in the sample image data set from the first plurality of features extracted from the sample image data set by means of a classification algorithm, wherein a second plurality of features extracted from at least one training image data set serves as a training data set of the classification algorithm.

Therein, the presence or absence of the pattern is known for each training image data set. In particular, the pattern refers to a structural change in the retina, from which the respective image data set has been acquired by OCT measurement. For example, the presence or absence of such a pattern may be determined independently from the method described herein by diagnostic methods, particularly fluorescence fundus angiography (FFA).

In particular, the value of a label is automatically determined a from the first plurality of features of the sample image data set by means of the classification algorithm, wherein the presence of the pattern in the sample image data set can be predicted by means of the value of the label.

In particular, the label may be a binary number, wherein a value of 1 indicates that the pattern is predicted in the sample image data set, and wherein a value of 0 indicates that the pattern is not predicted in the sample image data set.

In the context of the present specification, the term "classification" is used in its meaning known in the art of mathematics and statistics. It refers to the problem of identifying to which of a set of categories a new observation belongs, on the basis of a training set of data containing observations whose category membership is known. Referring to the present invention, an OCT image data set serves as an observation, wherein the presence of the pattern in the OCT image data set represents a first category and the absence of the pattern represents a second category.

In particular, the classification algorithm is used to determine whether the pattern is present in an OCT image data set. More particularly, the presence refers to a structural change in the retina, and the classification algorithm is used to determine whether the structural change has occurred in the retina. Even more particularly, the structural change in the retina corresponds to a successful retinal therapy outcome, most particularly a successful selective retinal therapy outcome.

In particular, the classification algorithm is selected from a decision tree learning algorithm, particularly a Random Forest algorithm, a support vector machine, a neural network, particularly a convolutional neural network, or a boosting algorithm. Therein, at least a subset of the extracted features is used to perform candidate splits in a decision tree learning algorithm.

In particular, the classification algorithm is a Random Forest algorithm.

In the context of the present specification, the term "Random Forest algorithm" refers to a decision tree learning algorithm, which uses bootstrap aggregating and a random subset method.

In the context of the present specification, the term "bootstrap aggregating" refers to the process of repeatedly selecting a random sample with replacement of the training set and fitting trees to these samples in a decision tree learning algorithm.

In the context of the present specification, the term "random subset method" refers to the process of selecting a random subset of features at each candidate split in a decision tree learning algorithm.

In certain embodiments, the method comprises directing a treatment laser beam at the retina. In particular, the treatment laser beam is adapted for retina therapy, more particularly selective retina therapy.

In certain embodiments, at least one treatment laser pulse of the treatment laser beam is applied to the retina.

In certain embodiments, the at least one treatment laser pulse comprises a pulse energy from 10 to 200 μύ, particularly from 40μύ to 140μύ.

In certain embodiments, a plurality of treatment laser pulses is applied to the retina.

In certain embodiments, the plurality of treatment laser pulses are applied with an inter- pulse interval of 2 to 20 ms, particularly 8 to 12 ms.

In certain embodiments, the treatment laser beam is coupled into the measuring arm of the OCT device.

In certain embodiments, the treatment laser beam is focused on the retina by means of an ophthalmoscope, particularly a digital scanning ophthalmoscope.

In certain embodiments, a statistical measure of the intensity values of at least a subset of the plurality of image data points from the optical coherence tomography image data set is extracted as a feature of the first plurality of features.

In certain embodiments, the statistical measure is a mean, a median, a root mean square, a variance, or a standard deviation. In certain embodiments, a mathematical function of the statistical measure of the intensity values of at least a subset of the plurality of image data points from the optical coherence tomography image data set is extracted as a feature of the first plurality of features.

In certain embodiments, the statistical measure is a mean, a median, a root mean square, a variance, or a standard deviation.

In certain embodiments, the intensity values of at least one subset of the plurality of image data points are aggregated into an aggregate value by means of an aggregation function, wherein a statistical measure or a function of the statistical measure of the aggregate value is extracted as a feature of the first plurality of features.

In certain embodiments, the aggregation function is a root mean square function or a sum function.

In certain embodiments, the intensity values of at least one subset of the plurality of image data points are aggregated into an aggregate value by means of an aggregation function, wherein a function of the statistical measure of the aggregate value is extracted as a feature of the first plurality of features.

In certain embodiments, the aggregation function is a root mean square function or a sum function.

In certain embodiments, the image data points of the subset of the plurality of image data points each comprise the same time coordinate value. That is, the intensity values of groups of image data points comprising the same time coordinate value are aggregated in to a respective aggregate value. Thereby, the two-dimensional sample image data set is transformed into a one-dimensional distribution, wherein the time coordinate value serves as an independent variable.

Such an aggregation of data allows a better analysis of time dependent effects and patterns in the sample image data set.

In certain embodiments, the image data points of the subset of image data points each comprise the same spatial coordinate value. That is, the intensity values of groups of image data points comprising the same spatial coordinate value are aggregated in to a respective aggregate value. Thereby, the two-dimensional sample image data set is transformed into a one-dimensional distribution, wherein the spatial coordinate value serves as an independent variable.

Such an aggregation of data allows a better analysis of axial-depth dependent effects and patterns in the sample image data set. In certain embodiments, the number of aggregate values exceeding a threshold subset value is extracted as a feature of the first plurality of features from the sample image data set.

In certain embodiments, the aggregation function is described by the formula

wherein S t is denoted as a group sum of intensity values X(l,t) calculated from groups of image data points having the same time coordinate value, I labels the spatial coordinate value, t labels the time coordinate value, and L labels the maximum number of spatial coordinate values of the sample image data set, wherein particularly a respective offset term is subtracted from the respective group sum, more particularly by a moving average operation, resulting in a plurality of corrected group sums, wherein the number of group sums or the number of corrected group sums exceeding a threshold value is extracted as a feature of the first plurality of features.

In certain embodiments, the sample image data set is split into a plurality of blocks of image data points, wherein the time coordinate values of the image data points in each block are in a respective partial range, and wherein a respective block feature of the first plurality of features is extracted from at least one block.

In the context of the present specification, the term block feature describes a feature, which is allocated to a respective block of image data points.

Dividing the sample image data set into blocks, allows a better analysis of periodic effects, particularly the effects of repeated pulses of a treatment laser beam, more particularly in retinal therapy, most particularly in selective retinal therapy.

In certain embodiments, the intensity values of at least one block are aggregated into a respective block aggregate value by means of an aggregation function, wherein a statistical measure of the block aggregate value is extracted as a block feature of the first plurality of features.

In certain embodiments, the aggregation function is a mean, a median, a root mean square, a variance, or a standard deviation.

In certain embodiments, a plurality of image data points from the sample image data set is allocated to a reference block, wherein for at least one image data point of each block a corrected intensity value is generated by subtracting from the respective intensity value of the at least one image data point of each block a reference intensity value of a corresponding image data point of the reference block. Subtracting the image data points of a reference block allows to analyze specific effects and patterns, which are not present in the reference block, particularly the effects of a treatment laser beam in retinal therapy, more particularly selective retina therapy.

In certain embodiments, at least one of the block features of the first plurality of features is extracted by means of the formula

wherein I labels the spatial coordinate value, t labels the time coordinate value, L n labels the maximum number of spatial coordinate values of the respective block, T n labels the maximum number of time coordinate values of the respective block, B ,n (l,t) is a reference-subtracted intensity value derived by subtraction of a respective intensity value taken from a reference block from the intensity value of the image data point in the block, and μ η ' is a mean determined by means of the formula

/½.' = r∑i t Tn B n ' (I, t). In certain embodiments, a maximum difference between subsequent block features of subsequent blocks is extracted as an additional feature of the first plurality of features. Therein the term subsequent blocks designates two blocks of image data points comprising subsequent time coordinate values.

In certain embodiments, at least one of the block features of the first plurality of features is extracted by means of the formula

wherein v n is a variance of the block B n determined according to the formula

and wherein v 0 is a variance of a reference block determined according to the formula

wherein I labels the spatial coordinate value, t and t' label the time coordinate values, L labels the maximum number of spatial coordinate values of the respective block, T n labels the maximum number of time coordinate values of the respective block, R n (t) is a root mean square of the intensity values of the image data points of the block B n at the time coordinate value t, and R 0 (t) is a root mean square of the intensity value of the image data points (x i 0 ) in the reference block, wherein R n and R 0 are respectively determined by means of the formula

In certain embodiments, a maximum value of the block features is extracted as an additional feature of the first plurality of features.

In certain embodiments, a frequency value, and/ or a phase value of at least one image data point is provided by means of a Fourier transformation, particularly a short time Fourier transformation, of the sample image data set, or by means of a spectrometer, particularly comprised in a system for predicting a pattern in an optical coherence tomography sample image data set according to the second aspect of the invention, wherein at least one feature of the first plurality of features is extracted from the frequency values and/ or the phase values.

A Fourier transformation allows to generate frequency values, and/ or phase values from the intensity values of image data points of an OCT sample image data set in the absence of a spectrometer.

In certain embodiments, a crest factor is determined from the intensity values, frequency values, or phase values.

The phase values in the OCT data contain information about sub-resolution

movements of scatterers as well as the velocity and movement pattern in which the displacement of scattering particles happens. This may provide information about the number, size and movement of microbubbles, the expansion and collapse of the RPE cell volume and/ or persisting sub-resolution changes in the layer structure.

In certain embodiments, a statistical measure, particularly a mean, a median, a root mean square, a variance, or a standard deviation, of the phase values is extracted as a feature of the first plurality of features.

In certain embodiments, a statistical measure of the differences of phase values between two time coordinate values is extracted as a feature of the first plurality of features.

This may indicate fast, random (chaotic) movement of scatterers and could hint at microbubbles, their number and size as well as their lifetime.

In certain embodiments, a feature of the first plurality of features is extracted from a time course of phase values at a given spatial position or a time course of an aggregate value, particularly a sum, an absolute sum, or a root mean square, of phase values, in particular wherein a maximum phase value or aggregate value, a difference between a maximum phase value or aggregate value and a minimum phase value or aggregate value, or a maximum difference between subsequent phase values or aggregate values is extracted as a feature of the first plurality of features. Therein, the term "subsequent phase values or aggregate values" designates pairs of phase values or aggregate values having a subsequent time coordinate value.

In certain embodiments, the number of aggregate values, particularly sums, exceeding 80 % of the maximum value is extracted as a feature of the first plurality of features.

In certain embodiments, a statistical measure, particularly a mean, a median, a root mean square, a variance, or a standard deviation, of the distribution along the time coordinate of aggregate values, particularly sums, of phase values is extracted as a feature of the first plurality of features.

Such features may indicate in what retinal layer effects occurred and to what extent the layers were affected.

In certain embodiments, a feature of the first plurality of features is extracted from a distribution of phase values having the same time coordinate value, wherein particularly a statistical measure, more particularly a mean, a median, a root mean square, a variance, or a standard deviation, is extracted as a feature.

In certain embodiments, a difference between the statistical measure, and a reference value is extracted as a feature of the first plurality of features.

In certain embodiments, a feature of the first plurality of features is extracted from a profile of phase values having the same time coordinate value along the spatial coordinate, particular wherein a maximum phase value, a difference between a maximum phase value and a minimum phase value, or a maximum difference between subsequent phase values is extracted as a feature of the first plurality of features.

Therein, the term "subsequent phase values" designates pairs of phase values having a subsequent spatial coordinate value.

In certain embodiments, a difference between the maximum phase value, difference between the maximum phase value and the minimal phase value, or the maximum difference between subsequent phase values, and a reference value is extracted as a feature of the first plurality of features.

This may indicate sub-resolution but persisting distortions in the scatterer distribution. In certain embodiments, a plurality of differences between subsequent phase values along the time coordinate are calculated, and a sum of the differences is extracted as a feature of the first plurality of features.

In certain embodiments, a plurality of differences between subsequent phase values along the spatial coordinate are calculated, and a sum of the differences is extracted as a feature of the first plurality of features.

In certain embodiments, a plurality of first differences between subsequent phase values along the time coordinate are calculated, and a plurality of second differences between subsequent phase values along the spatial coordinate are calculated, and the sum of the first differences and the second differences is extracted as a feature of the first plurality of features.

In certain embodiments, a plurality of differences between subsequent phase values along the time coordinate are calculated, and a statistical measure, particularly a mean, a median, a root mean square, a variance, or a standard deviation, of the differences is extracted as a feature of the first plurality of features.

In certain embodiments, a plurality of differences between subsequent phase values along the spatial coordinate are calculated, and a statistical measure, particularly a mean, a median, a root mean square, a variance, or a standard deviation, of the differences is extracted as a feature of the first plurality of features.

In certain embodiments, a plurality of first differences between subsequent phase values along the time coordinate are calculated, and a plurality of second differences between subsequent phase values along the spatial coordinate are calculated, and a statistical measure, particularly a mean, a median, a root mean square, a variance, or a standard deviation, of the first differences and the second differences is extracted as a feature of the first plurality of features.

In certain embodiments, a spectrogram comprising a plurality of time and frequency dependent amplitude values is generated by means of a Fourier Transformation, particularly a Short-Time-Fourier-Transformation (STFT) according to the formula

wherein ω is a frequency, t labels a first time coordinate of the image data points, h(t - τ) is a windowing function, T labels a second time coordinate of the windowing function, R(x) is an aggregate value generated by means of a root mean square aggregation function. Therein, the spectrogram is converted to a sum spectrogram ( ) by means of the formula Ζ(τ =∑ ω Ζ(τ, ω ,

wherein particularly the amplitude values of the spectrogram are median filtered prior to conversion to the sum spectrogram, and the sum spectrogram is normalized, particularly by means of the formula Ζ(τ) = Z(T) - min(Z(T)) ,

max (Z ( -min (Z (

wherein max (Ζ(τ)) is the maximum of the sum spectrogram, and min(Z(x)) is the minimum of the sum spectrogram, resulting in a normalized spectrogram. Therein, a crest factor calculated according to the formula

max(Z(r))

u =—____==

^Ψ∑ Τ τ=ί \ ζ ω\ 2 is extracted as a feature of the first plurality of features, wherein T labels the maximum number of time coordinate values of the normalized spectrogram ζ(τ), wherein particularly the amplitude values of the spectrogram are converted to a dB scale prior to generating the crest factor.

In certain embodiments, the sum ζ(τ) is determined as an additional feature of the first plurality of features.

According to a second aspect of the invention, a system for predicting a pattern in an optical coherence tomography sample image data set is provided. The system comprises an optical coherence tomography (OCT) device comprising a first laser source, a measuring arm, and a reference arm. The first laser source, the measuring arm, and the reference arm are optically connected, particularly by means of a fiber coupler.

The system further comprises an ophthalmoscope, particularly a scanning digital ophthalmoscope (SDO), which is coupled into the light path of the measuring arm of the optical coherence tomography (OCT) device, such that a plurality of light intensities of reflected light, and/or diffracted light, and/or a plurality of phases in a retina of an eye observed by the ophthalmoscope are measurable by means of the optical coherence tomography device.

The system further comprises a second laser source, which is coupled into the light path of the measuring arm of the optical coherence tomography device, particularly by means of a dichroic mirror, such that a retina of an eye is observable by means of the optical coherence tomography device during illumination of the retina by the second laser source. The system is characterized by an analysis device, which is adapted to automatically extract a first plurality of features from the sample image data set received from the optical coherence tomography (OCT) device, and automatically predict the presence or absence of a pattern in the sample image data set from the first plurality of features by means of a classification algorithm, wherein the pattern corresponds to a structural change in a retina of an eye, and wherein a second plurality of features extracted from at least one training image data set serves as a training data set of the classification algorithm, particularly by means of steps of the method according to the first aspect of the invention.

In certain embodiments, the system comprises a spectrometer, which is optically connected to the measuring arm, particularly by means of the fiber coupler, such that a frequency spectrum of the reflected and/ or diffracted light is measurable by means of the spectrometer.

According to a third aspect of the invention, a computer program is provided. The computer program comprises program code, which is adapted to automatically extract a first plurality of features from a sample image data set received from an optical coherence tomography (OCT) device, and automatically predicting the presence or absence of a pattern in the sample image data set from the first plurality of features by means of a classification algorithm, wherein the pattern corresponds to a structural change in a retina of an eye, and wherein a second plurality of features extracted from at least one training image data set serves as a training data set of the classification algorithm, when the program is executed on a computer. In particular, the program code is adapted to perform steps of the method according to the first aspect of the invention when executed on a computer.

The invention is further illustrated by the following examples and figures, from which further embodiments and advantages can be drawn. These examples are meant to illustrate the invention but not to limit its scope.

Short description of the figures

Fig. 1 shows a system according to the invention comprising an optical

coherence tomography (OCT) device, a second laser source, and an analysis device.

Fig. 2 shows a time sequence of the application of five laser pulses and the corresponding OCT M-scan data. Fig. 3 shows an overview of the proposed features, including blockwise analysis of M-scans, speckle variance features and spectrogram features.

Fig. 4 shows Color Fundus and Fundus Fluorescein Angiography (FFA)

images of an eye region after selective retina therapy,

Fig. 5 shows 60 ms extracts of examples for manually labeled M-scans.

Fig. 6 shows receiver operating characteristic (ROC) curves generated by means of the algorithm of the present invention on ex-vivo data,

Fig. 7 shows the performance of the algorithm according to the invention

depending on the duration of laser application and the subdivision of an

M-scan.

Fig. 8 shows a performance analysis of the classification algorithm on clinical in-vivo Selective retina therapy (SRT) data. Detailed description of the figures

Fig. 1 shows a system 1 for predicting a structure S in an optical coherence

tomography sample image data set X, wherein the system 1 comprises an optical coherence tomography (OCT) device 2, a second laser source 3, particularly a laser source adapted for retina therapy, more particularly selective retina therapy, and an analysis device 4.

The depicted OCT device 2 comprises a first laser source 21 , a measuring arm 22, and a reference arm 23, wherein the first laser source 21 , the measuring arm 22, and the reference arm 23 are connected by means of a fiber coupler 25, and fiberoptic connections F, such that light emanating from the first laser source 21 is directed to the measuring arm 22, and to the reference arm 23 by means of the fiberoptic connections F and the fiber coupler 25.

Therein, the measuring arm 22 is defined as the light path between the fiber coupler 25 and the retina 51 of the eye 5, and the reference arm 23 is defined as the light path between the fiber coupler 25 and the reflective element 231 . In particular, the reflective element 231 may be a mirror.

An alternative OCT device 2 according to the invention may be designed such that light emanating from the first laser source 21 is directed to the measuring arm 22, and/or to the reference 23 arm without any fiberoptic connections F, particularly by means of at least one beam splitter and/or at least one mirror. According to a first light path, light emitted by the first laser source 21 travels through the fiber coupler 25 and the measuring arm 22, is reflected and/ or diffracted by the retina 41 , travels back to the fiber coupler 25, and to the detector 242 of the

spectrometer 24.

According to a second light path, light emitted by the first laser source 21 travels through the fiber coupler 25 and the reference arm 23, is reflected by the reflective element 231 , is condensed by the lens 232, travels back to the fiber coupler 25, and to the detector 242 of the spectrometer 24.

The phase difference of light of the first light path and the second light path is determined by the difference of the length of the measuring arm 22 and the reference arm 23. If the measuring arm 22 and the reference arm 23 have the same length, a maximal constructive interference between the light beams of the measuring arm 22 and the reference arm 23 occurs. In particular, the length of the reference arm 23 can be changed by altering the position of the reflecting element 231. The total light intensity as a result of interference between the first light path and the second light path is measured as an OCT signal by the detector 242.

In particular, the length of the measuring arm 22 depends on the position in the axial direction a, and the refractive index of the section of the retina 51 illuminated by the light beam. Therefore, topological and structural changes in the retina can be measured by means of the OCT signal.

An ophthalmoscope 27, particularly a digital scanning ophthalmoscope is integrated into the light path L of the measuring arm. By means of the ophthalmoscope 27, the light beam is focused in an axial direction a, particularly such that a specific layer of the retina 51 is in the focus of the light beam.

Fig. 1 further depicts a scanner 28, which is adapted to change the direction of the light beam particularly in a horizontal direction h, and/or in a vertical direction v, such that a specific region of the retina 51 is illuminated by the light beam. In particular, the scanner 28 is adapted to change the direction of the light beam in both the horizontal direction h and the vertical direction v.

In the spectrometer 24, the light beam is split up into light of different wavelengths by diffraction on the diffractive element 241 , and light of the different wavelengths is detected by the detector 242. In particular, the diffractive element 241 may be a grating. In particular, the detector may be a charge coupled device (CCD) detector.

According to an alternative embodiment, no spectrometer 24 is integrated into the OCT device 2. Therein, light of the measuring arm 22, and the reference arm 23 particularly travels directly from the fiber coupler 25 to the detector 242 without being split up into different wavelengths.

Fig. 1 further shows a second laser source 3, wherein the light beam emitted by the second laser source 3 is directed to a dichroic mirror 28 positioned in the measuring arm 22 between the scanner 26 and the ophthalmoscope 27. The light beam emitted by the first laser source 21 and the light beam emitted by the second laser source 3 can be combined by means of the dichroic mirror 28. Therein, the light beam emitted by the second laser source 3 can be coupled into the light path of the measuring arm 22.

The light beam emitted by the second laser source 3 can thus be directed to the ophthalmoscope 27 and particularly focused onto the same part of the retina 51 that is illuminated by the first laser source 21 .

In particular, the second laser source 3 is adapted to emit a light beam of a power sufficient for selective retina therapy. Therein, light pulses emitted by the second laser source 3 particularly cause evaporation of water in cells of the retinal pigment epithelium, resulting in the formation of microbubbles.

Furthermore, an analysis device 4 is shown in Fig. 1. The analysis device 4 is connected to the detector 242, such that image data may be extracted from the detector by means of the analysis device 4.

In particular, the analysis device 4 may comprise a microcomputer, or an interface adapted to extract image data from the detector 242 and transmit the image data to a computer.

In particular, software running on the analysis device 4 is adapted to execute steps of the method according to the invention.

Fig. 2 (a) shows a time sequence of the application of five treatment laser pulses over a time-course 60 ms. Therein, a peak indicates that the treatment laser was switched on during the respective time interval. The depicted laser pulses emanate from a second laser source 3 of a system depicted in Fig. 1 or an analogous system, and are applied to a retina 51 of an eye 5, particularly during selective retina therapy (SRT) by means of the system depicted in Fig. 1 or an analogous system. Fig. 2 (b) shows the corresponding optical coherence tomography (OCT) M-scan data. Therein, the term IVI- scan refers to a time series of A-scans, and the term A-scan refers to an optical coherence tomography scan, wherein the sample is scanned in an axial direction a. During a scan in an axial direction a (A-scan), the focus of the used optical system, particularly an ophthalmoscope, is changed such that during each OCT measurement the OCT laser beam, particularly the laser beam emanating from the first laser source 21 shown in Fig. 1 , is reflected and/ or diffracted by a specific layer of the retina.

The y-coordinate of the plot depicted in Fig. 2 (b) represents the axial direction a, wherein the arrow designates the direction of the scan from the surface of the retina to the photoreceptor cells. A scale bar in the lower left corner depicts the length of the y- The x-coordinate represents the time course of A-scan measurements constituting the M-scan. The grey values of the pixels of the plot represent the intensity of the OCT signal, wherein light grey values represent higher intensity and dark grey values represent lower intensity. In particular, high OCT signal intensities are a result of constructive interference of the light beams in the sample arm and the reference arm, and depend on the distance of the respective illuminated specimen from the OCT system, the intensity of light reflected/ diffracted from the sample, the corresponding phase of the reflected/ diffracted light, and the length of the reference arm.

The label "PRL" indicates the position of the pre-retinal layer in the retina, the label "RPE" indicates the position of the retinal pigment epithelium, which represents the layer which is targeted by selective retina therapy (SRT), and the label "CHR" indicates the position of the choroidea.

The white arrows at the top of the plot indicate the time points at which the treatment laser was switched on. These time points correspond to the peaks depicted in Fig. 2 (a).

Fig. 3 shows a schematic representation of the feature extraction process from an OCT M-scan data set X. The OCT M-scan data set X consists of a plurality of image data points x, which constitute a speckle pattern SP. In particular, features may be extracted from the OCT M-scan data set X by four different methods. In blockwise M-scan feature extraction F1 , the OCT M-scan sample image data set X is split into a plurality of blocks B n , and features are extracted from each block B n . Likewise, in blockwise speckle feature extraction F2, the OCT M-scan sample image data set X is split into a plurality of blocks B n . However in the case of blockwise speckle feature extraction F2, the root mean square (rms) is first computed individually for all temporal sampling points and the root mean square values are analyzed for each block. In the extraction of speckle variance features F3, the root mean square (rms) is first computed individually for all temporal sampling points and the root mean square values are analyzed for the entire OCT M-scan data set X. In the extraction of spectrogram features F4, a spectrogram is extracted from the time series of the OCT M-scan data set X, particularly by means of a Short-Time Fourier Transform (STFT). Fig. 4 shows SRT lesions in (a) Color Fundus image and (b) Fundus Fluorescein Angiography (FFA) for the same eye region. With appropriate laser energies, lesions after selective retina therapy remain invisible in Color Fundus image while being visible in FFA (state-of-the-art for identification of lesions introduced by SRT.

Fig. 5 shows 60 ms extracts of examples for manually labeled M-scans. Fig. 5(a) shows an ex-vivo OCT M-scan with no detectable signal variations labeled as class "0". Fig. 5 (b) and (c) show M-scans with detectable signal variations limited to the RPE / Bruch's membrane complex and throughout the retina, respectively. Scans (b) and (c) were consequently labeled as class "1 ".

Fig. 6 (a) shows an ROC curve of our algorithm on ex-vivo data. The AUC is 0.996 and the arrow indicates the threshold value for 95% specificity used for clinical evaluation of the classification performance. Fig. 6 (b) shows ROC curves showing the performance comparison for different subsets of used features (ex-vivo data).

Fig. 7 shows algorithm performance depending on the duration of laser application (a) and the subdivision of an M-scan (b).

Fig. 8 shows a performance analysis of the classification algorithm on clinical in-vivo SRT data. Classifications were evaluated for 100% and 95% specificity as well as for 100% sensitivity. Classification results and the manual labeling were compared to the FFA visibility as assessed by the attending ophthalmologist. Classification performance was best for 95% specificity, false results are highlighted.

Example

Purpose In recent years, selective retina laser treatment (SRT), a subthreshold therapy method, avoids widespread damage to all retinal layers by targeting only a few. While these methods facilitate faster healing, their lack of visual feedback during treatment represents a considerable shortcoming as induced lesions remain invisible with conventional imaging and make clinical use challenging. To overcome this, we present a new strategy to provide location-specific and contact-free automatic feedback of SRT laser applications.

Methods We leverage time-resolved optical coherence tomography (OCT) to provide informative feedback to clinicians on outcomes of location-specific treatment. By coupling an OCT system to SRT treatment laser, we visualize structural changes in the retinal layers as they occur via time-resolved depth images. We then propose a novel strategy for automatic assessment of such time-resolved OCT images. To achieve this, we introduce novel image features for this task, that when combined with standard machine learning classifiers yield excellent treatment outcome classification capabilities.

Results Our approach was evaluated on both ex-vivo porcine eyes, as well as human patients in a clinical setting, yielding performances above 95% accuracy for predicting patient treatment outcomes. In addition, we show that accurate outcomes for human patients can be estimated even when our method is trained using only ex-vivo porcine data.

Conclusions The proposed technique presents a much needed strategy towards noninvasive, safe, reliable and repeatable SRT applications. These results are

encouraging for the broader use of new treatment options for neovascularization based retinal pathologies.

Selective retina therapy (SRT) is an efficient laser treatment for a variety of retinal pathologies. At its core, the laser inducing a photo-mechanic disruption by evaporating water in targeted cells. In contrast to conventional laser photocoagulation, which is prone to induce laser damage to all retinal layers, SRT achieves therapeutic effects by solely targeting the retinal pigment epithelium (RPE) which is responsible for nourishing light absorption cells. By doing so, SRT not only reduces negative side effects of traditional photocoagulation but also facilitates healing of neighboring tissue. In particular, recent studies have shown promising long term efficacy of SRT in both animal and patient studies when the treatment is executed with laser energies that produce visible lesions in angiography, but invisible to ophthalmoscopic imaging, as illustrated in Fig. 4.

While highly beneficial, SRT suffers from the main drawback of missing direct visual feedback regarding the success of the therapy. Introduced lesions remain invisible in conventional ophthalmoscopic imaging due the absence of coagulation in the inner retinal layers (see Fig. 4). For this reason, the selection of an appropriate laser energy and a reliable monitoring of the therapy are both challenging and critical. The patient- specific concentration of melanin in the RPE further aggravates the determination of the necessary treatment energy level as it influences the rate of conversion from laser energy to heat. For these reasons, real-time and objective evaluation of the introduced tissue damage as it occurs is paramount for safe, reliable and repeatable SRT.

Current approaches to SRT monitoring are either based on the invasive and time- consuming fundus fluorescein angiography (FFA) or the detection of ultrasonic pressure waves of collapsed cells or change analysis in the backreflection during the presence of laser-induced microbubbles. While the latter methods have already been implemented, these approaches suffer from the lack of depth information. Recent research has shown however that optical coherence tomography (OCT), acquired simultaneously with the laser therapy may provide the missing spatial and temporal information necessary for a comprehensive, repeatable and reliable lesion assessment. In a pilot study, the application of single SRT laser pulses to induce RPE lesions was detectable in OCT data and appeared to correlate well with the extent of tissue damage imaged with FFA.

This visible FFA leakage is a consequence of the accumulation of contrast agent in sub-retinal tissue which breached the blood retina barrier. Yet no method to date has been able to automatically assess laser treatment outcome at specifically targeted locations. This inability severely hinders the clinicians' capacity to use SRT as a treatment option.

To overcome this limitation, this example introduces the first automatic and observer- independent classification scheme for time-resolved OCT data of SRT lesions. To this end, we introduce novel image features for time-resolved OCT that when combined with traditional classification schemes, provide excellent identification of positive and negative treatment outcomes. The proposed approach was evaluated on both ex-vivo porcine eye samples and human patients undergoing SRT. In addition, we

demonstrate here that our features are reliable for classification, to the point that classification models trained on ex-vivo porcine data can effectively be used for prediction in human patients. As such we show here that use of time-resolved OCT during SRT can provide a direct feedback to the ophthalmologist and allows SRT to be executed using strongly reduced pulse energies

without the risk of undertreatment.

System Overview In this section, we describe our SRT-OCT system used for patient treatment. The OCT system used for data acquisition shown in Fig 1 features a line scan frequency of 70 kHz and a spectral bandwidth of 170 nm centered at 830 nm (EBS8C10, Exalos AG, Switzerland), leading to an axial resolution of 1.78 μηη in air and is described in detail elsewhere (Steiner et al., Retinal laser lesion visibility in simultaneous ultra-high axial resolution optical coherence tomography. IEEE Photonics Journal 6(6), 1-1 1 (2014)).

While a Fourier-domain OCT system was used for the experiments in this paper, any type of OCT system (e.g. time-domain, swept source) could be used instead, as long as the system parameters are properly set. The OCT was combined with the treatment system using a dichroic mirror to enable simultaneous data acquisition and M-scans (i.e. repeated A-Scans at the same XY position over time) were acquired during each SRT pulse application.

Prior to any further analysis and processing, the raw OCT data undergoes conventional pre-processing steps such as λ-k mapping, numerical dispersion correction, background subtraction, mapping to a color map and axial motion correction.

The SRT pulses are executed using a frequency-doubled Nd:YLF laser (SRT Vario, Medical Laser Center Lubeck, Lubeck, Germany) with pulse width of 250 ns and a pulse repetition rate of 100 Hz for 30 pulse-trains (see Fig. 2(a)).

Laser energy was guided onto the retina using a scanning digital ophthalmoscope (SDO, Wild medTec, Austria) and SRT lesions were placed onto the retina by the ophthalmologist using a built-in manual detection mirror. OCT M-scans are acquired as a single depth profile at a fixed position on the sample. As such, M-scans provide a time-resolved depth profile of light scattering and reflection in the tissue under investigation (see Fig. 2(b)). The temporal resolution thereby depends on the maximal A-Scan rate of the OCT system and the axial resolution is defined by the coherence length of the employed OCT system. With this, M-scans are capable of mapping even sub-resolution changes in the scattering structure by abstract variations in the phase of the OCT signal. For SRT, where the laser pulse parameters circumvent thermal effects, the damage in the RPE cells is based on photomechanical rather than photothermal effects.

Thus, detectable features in OCT M-scans may primarily be based on the presence of microbubbles in the RPE cells due to the induced water evaporation. Extensive recent studies have confirmed the presence of a set of distinctive signal variations in the recorded OCT data that correlate well with lesion visibility in fundus fluorescein angiography. It is thus assumed that the presence of temporal OCT signal variations as shown in Fig. 2(b) represents successful RPE cell rupture.

Problem Statement In this example, we consider the problem of automatically identifying the successful (i.e. the presence of RPE cell rupture) and unsuccessful SRT laser application from the corresponding M-scan image data. Let X k , X k , k = 1,2, ... , K, be an L x T M-scan- OCT image of the kt observed SRT laser application with L and T representing the spatial and temporal dimensions respectively. Then, let y fe e {0,1) be the label associated with the corresponding X k , where y k = 1 if the therapy is successful, and where y k = 0 otherwise

We define a feature extractor function g-. R LxT → R d which maps an l x T image X k to a d-dimensional feature vector. Our aim then is to learn a classifier /: R d → {0,1) assigning a binary label to a given d-dimensional vector. In our case, a Random Forest classifier is used to represent our function f. We now introduce a novel feature set g to tackle this classification problem.

M-scan Features In this Section, we describe our approach to the classification problem presented above. Accordingly, we describe how representative features can be extracted from M-scans. Successful SRT application becomes visible in M-scans as high frequency intensity variations (Fig. 2(b)), in standard OCT images during and following the single SRT pulse application.

The physical principle of OCT imaging only allows the detection of the reflectivity and the position of reflective components in the sample tissue. As such, detectable features are limited to variations in intensity, speckle pattern and the phase of the OCT signal.

Features were thus computed from pixel intensity-based speckle images, blockwise image analysis as well as from time-frequency analysis (TFA) using a Short-Time Fourier Transform (STFT). Fig. 3 shows an overview of the features we compute here and which we now discuss in detail.

Blockwise M-scan Features To represent variations in the temporal pixel intensity distribution of the OCT M-scan, a blockwise speckle analysis is performed by dividing the OCT M-scan into equal size-blocks of signal. The blocks are defined such that only one of the two subsequent blocks will contain the temporal position of laser application as shown in Fig 3.

Let us denote the nth block as B e R Lt > xTb , n = 1, ... , N— 1, and the reference block B Q e ff L6XT6 which are taken from the kt M-scan X K . Here, T B and L B show temporal and spatial dimensions of each block respectively. Note that B corresponds to image data before any laser application has been performed.

Henceforth, we drop the subscript/superscript k from the notation for the sake of simplicity. As the signal variations inside a defined block can be inferred from the variance of the pixel distributions, the follo — 1) dimensional feature vector is computed as in the following:

where B n ' = B n - B 0 is the reference-subtracted block and μ ' = b t ' b B n ' l, t)is its mean. Moreover, the maximum gradient of the resulting vector u m = max(Vu i ) is also used as a feature. The blockwise M-scan feature vector then contains the two extracted feature components, i.e. u bm t = [ui m t , u m ] of size N.

Blockwise Speckle Features While large movements alter the spatial intensity distribution in the OCT images, smaller effects may only be detectable in the speckle pattern which, in contrast to the common intensity image representation, is sensitive to sub-wavelength vibrations. In order to provide a broad set of features for classification, time resolved variations in speckle pattern were calculated from the OCT M-scans. To further emphasize signal blocks containing high signal variations, the root mean square (rms) is first computed individually for all temporal sampling points leading to a vector with length T. Subsequently, the variance of the vector values within each signal block is computed for all N blocks, i.e.,

V n = ) | 2 IS the rms of nth block at temporal position t and v n corresponds to the variance of such rms values within the nth block computed along time axis. Then, an (N— 1) -dimensional blockwise speckle feature vector is extracted by normalizing the computed rms variance vector as in the following: ... , N— 1, where v 0 holds for the rms variance of the reference block. In addition, we include the maximum value of the resulting vector, that is u s (n) = max(u^ s ). Accordingly, one gets a blockwise speckle feature vector u bs t = [ u fe t j U 6s] of S j ze ^v.

Speckle Variance Features A second set of features based on OCT speckle information is computed based on the temporal variation of the speckle values in the OCT signal. For this, we identify the variations along time by the following sum:

S t Χ(1, t), where S(t) represents the intensity sum for each time step t for the OCT image X. In order to eliminate the bias due to the varying background, the offset term S 0 //(t) for the corresponding observation is subtracted in the following:

S(t)' = S(t) - S 0 ff (i), in which S(t) and S 0 //(t) represent zero-mean and unit variance correspondences of S(t) and S of f (i) respectively. The offset term above is computed by the moving average operation. In the resulting vector S(t)' , time-resolved variations can be quantified with the number of values surpassing an empirical threshold s th as shown by:

usv _ £T = i i{- seJj | S>s } (S(t)') > where I is the indicator function and u sv is the resulting feature related to the speckle variance.

Spectrogram Features It is believed that the presence of microbubbles in the RPE layer leads to rapid rearrangement of the scattering structure of the retinal tissue which is represented as high-frequency variations in the phase and intensity distribution of the OCT signal. We therefore introduce a time-frequency analysis by means of a Short- Time Fourier-Transform (STFT) to encode such signal variations in the feature set.

Here, STFT is performed on the rms signal ff (£) of an OCT sample X as

Ζ(τ, ω) =∑ =1 R( )h(t— τ)θχρ(-/ωτ) , with h( ) and (ω, τ) being the windowing function and the STFT of rms signal respectively. We then perform a median filtering operation on the amplitude of the resulting STFT and call \t Z m ( , ω). Subsequently, the sum of the median filtered spectrogram magnitudes is found for each time component, i.e.,

Z(T) = Σω(. ζ πι( τ > ω ))άΒ < where (. ) dB subscript represents the conversion to the dB scale. With the normalization given by:

Z( ) =— ζ ( τ _)-" ι η ( ζ (τ)) we com p Ute crest f ac tor in the following

J max(Z(i)-mm(Z(i)) ^ 3

= num b er of the evaluated signal windows by STFT

and which constitutes the first component of the spectrogram feature. Finally, the global description of the spectrogram is encoded with the sum u 2 sp =∑ =1 Ζ(τ), which is also added to the feature vector. All spectrogram features are then combined to create the two-dimensional spectrogram feature vector u sp t = [u^ u^] .

As a final step, blockwise M-scan, speckle variance and the spectrogram features are concatenated in one vector, i.e. = [u bm t , u bs t , u sv , u sp t ], which constitutes our final feature vector of length 2N + 3.

Experimental Results

We now detail the experimental evaluation of our automatic approach compared to manually evaluated M-scans. We show results on both ex-vivo and human clinical data.

Experimental Setup Two different datasets, namely ex-vivo and in-vivo human patient data were collected for the purpose of validation. Enucleated ex-vivo porcine eyes were collected from a local slaughterhouse and stored in DMEM solution. M-scan OCT images were recorded from 153 lesions from 22 enucleated porcine eyes with SRT pulse energies ranging from 1 10 μ J to 200 and a retinal treatment spot size of 170 nm in order to best represent effects from under- to over-treatment. Scans were acquired at the full speed of 70 kHz to provide the highest possible temporal resolution with incident OCT laser power of 0.7 mW. In addition to the ex-vivo data, 16 in-vivo human samples were recorded from 4 patients undergoing clinical SRT. Ethics approval (KEK-Nr. 003/10) to conduct this study was obtained from the local ethics committee, which is working in accordance with ICH-GCP guidelines. Patients were chosen and informed about the study according to the approved study protocol and patient consent was obtained

beforehand. For the treatment, applied pulse energies varied from 40μύ to 140μύ.

Before classification, all OCT M-scans underwent standard OCT post-processing as described above. To correct for axial motion, the lag of each A-Scan with regard to a reference A-Scan was determined by calculation of the cross-correlation and A-Scans were shifted accordingly. The application of the SRT laser was detected using a photodetector sampled at 1 MHz which allowed the determination of the temporal position of laser applications within the M-scan OCT data stream.

Using the information of the photodetector, M-scans were cropped to the length of the laser pulse train application of 300 ms. This resulted in M-scans of 600 x 20000 pixels. The ex-vivo M-scans were manually annotated (i.e. successful or unsuccessful treatment) by two independent observers based on the presence of high-frequency signal variations in the M-scan OCT images as depicted in Fig. 5. Unfortunately, in the case of full eye ex-vivo samples, no alternative evaluation method (e.g. histological evaluation) is available, as they suffer from artifacts or require significant preparation of the samples.

M-scans with inconsistent labels (9 out of 162) were excluded from further

consideration. From this, the training set consisted of a total of 153 M-scans from 22 ex-vivo eye samples including 1 19 positive and 34 negative M-scans. The human in- vivo dataset consisted of 16 clinical samples from 4 different patients. Unlike ex-vivo samples, in-vivo samples were labeled by an attending ophthalmologist using FFA, from which an assessment regarding each laser application was made and served as groundtruth.

Classifcation Performance As described above, block-wise, speckle and

spectrogram features were first extracted from M-scans. The length of block was set to 353 pixels with a pulse train of 100 Hz and the data acquisition at 70 kHz. The height of the blocks was 300 pixels centered at the RPE. With 30 pulses applied per treatment, this resulted in a total of 60 blocks. For spectrogram features, a Hamming window of length 189 and a hop factor of 9 was experimentally defined. Accordingly, a vector of 123 features was built and a Random Forest classifier of 200 trees was built and cross- validated with 50 iterations. In the following, we present the classification performance for ex-vivo and in-vivo samples.

Ex-vivo performance The training and testing of the classifier were performed first on ex-vivo samples. 102 OCT images out of 153 were used for training and the rest was used for testing. Fig. 6(a) depicts the ROC curve of our approach, with an Area-Under- Curve (AUC) of 0.996. Our results show that a 0.05 false positive rate (e.g. 1 in 20) still yields over a 0.95 true positive rate (e.g. 19 in 20) when detecting success treatment applications.

Feature Contribution In order to assess the impact and importance of the different features used by our classifier, classifiers were retrained using either only spectrogram features 61 or only speckle variation features 62, as well as all features combined 63. All other parameters such as numbers of trees, depth and number of iterations were kept constant. Fig. 6(b) shows the classification performance for these different settings.

As can be seen from these results, the speckle variation features 62 contribute the most to the classification. The spectrogram features 61 alone would already enable reasonable performance but it is only in combination with the others 63 that a high sensitivity and specificity can be attained. This performance may thus be further increased by including phase information and polarization sensitive features.

Temporal Features Distribution One of the ultimate goals of an automatic

classification treatment application is the possibility to cease the laser application once the desired effect is attained. In this sense, the temporal distribution of the features is relevant to stopping treatment early. To test the relevance of each part of the M-scan signal, the algorithm was trained with the features computed identically but with reduced signal length. The M-scans of 300 ms laser application were divided into three sets of 100 ms ("start" 71 , "middle" 72, and "end" 73), and algorithm performance was evaluated for each single set.

Fig. 7 shows the results of this analysis and presents the ROC for the each single 100 ms set as well as for the first 200 ms, consisting of the "start" and "middle" set 74, and the full set 75 of the entire 300 ms.

As can be seen in Fig. 7, the described features are much more prominent within the first two thirds of the applied pulse train. This may be caused by the break-up of melanin clusters occurring only at the beginning of the pulse application while later pulses merely lead to a rearrangement of the melanin complexes. As the temporal distribution of the features is highly interesting for the extension of the algorithm to a point where it may predict therapy outcome allowing to cease the laser once laser application was successful, future research must focus on the link between the quantitative amount of RPE damage and the extent of OCT signal variations.

ln-vivo Performance In a final step, we trained our algorithm using all the ex-vivo M- scans and keeping the same parameters as in the ex-vivo case; then tested our classifier on the in-vivo samples. The classification results were compared to both FFA images and manual OCT labeling. Fig. 8 shows the results on the 16 clinical SRT applications with FFA leakage evaluated by the attending ophthalmologist. The algorithm performance was evaluated at threshold values for 100% specificity, 95% specificity and 100% sensitivity. The best performance was thereby achieved using threshold for 95% specificity leading to a correct prediction of SRT outcome in 15 out of 16 cases corresponding to a success rate of 93.8% when compared to the FFA analysis.

With the classification tuned to prevent false positive results, the algorithm classified one lesion as positive where no leakage was detected in the FFA. However, visual inspection of the OCT scan revealed the presence of the described features, it thus remains unclear, whether the result is really a false-positive classification or not. For this corresponding lesion, pathological features in the area of treatment may have covered or omitted visible leakage in the FFA. There is, unfortunately, no way to further determine the true outcome of the treatment. The result was thus kept as a false classification.

Despite the fact that the algorithm in this paper was trained exclusively on ex vivo data, its performance on in-vivo human data is remarkably high. This can be somewhat explained by the fact that SRT omits heat dissipation and convection at too short time scales so that cooling effects due to blood perfusion do not play any role in the assessment. Nevertheless, this is a significant outcome as further studies may be based on enucleated porcine eyes which include no ethical considerations and are widely available. Moreover, we observed that the algorithm correctly classified also one lesion as a positive RPE rupture which was in agreement with the FFA while the visual inspection of the OCT data revealed no features. This finding may be a sign of the fact that the algorithm shows a higher specificity than visual inspection by human annotators.

Conclusion In this example, we have introduced a novel strategy for automatic classification of SRT lesions using M-scan OCT data. The presented features and algorithm showed encouraging performance and provided classification results in good agreement with clinical evaluation of lesions by FFA. In particular, our approach provided correct prediction of RPE leakage in 15 out of 16 cases in a preliminary human clinical trial. The proposed OCT-SRT system associated with a classification approach showed good performance for ex-vivo samples and we demonstrated that our approach can be trained using ex-vivo data even when evaluating in-vivo data.

In addition, the performance of the algorithm was found to be independent of the lateral position of the measurement spot within the treatment area as long as the

measurement and treatment area are fully overlapped.

Such OCT-based feedback and lesion classification may be extended to a predictive classification which can be used to cease the laser once the desired effect is achieved and may serve as an additional safety net for sub-threshold therapies.

Moreover, with the inclusion of volumetric data, captured before and after laser application, the algorithm performance may further be improved, enabling a more detailed classification. The use of OCT data has the potential of providing treatment feedback which is crucial for a reliable and repeatable therapy and would eventually support the proliferation of the reliable and cost-effective SRT as a treatment option.

List of reference signs

System for predicting a structure in a sample image 1 data set

Optical coherence tomography (OCT) device 2

First laser source 21

Measuring arm 22

Reference arm 23

Reflective element 231

Lens 232

Spectrometer 24

Diffractive element 241

Detector 242

Fiber coupler 25

Scanner 26

Ophthalmoscope 27

Dichroic mirror 28

Second laser source 3

Analysis device 4

Light path L

Fiberoptic connection F

Eye 5

Retina 51

Speckle variation features 61 Spectrogram features 62

Combined features 63

Start section 71

Middle section 72

End section 73

Start + middle section 74

Full scan 75

Structure S

Speckle pattern SP

Sample image data set X

Axial direction a

Horizontal direction h

Vertical direction v

Blockwise M-scan feature extraction F1

Blockwise speckle feature extraction F2

Speckle variance feature extraction F3

Spectrogram feature extraction F4