• Ei tuloksia

Publication III: Approximate marginalisation of absorption

OPTICAL TOMOGRAPHY

The fDOT image reconstruction is most often carried out using the so-called normalised Born approximation model [26], where the measure-ment vector is the measured fluorescent emission data vector yfobs scaled by the measured excitation data yeobs in (2.21). From the practical point of view, a significant feature of the Born normalised model is that it can tolerate inaccurately known target absorption and scattering distributions (µa,µs)to some extent. In practical experiments, the actual values of these parameters are usually not known accurately, and therefore one is bound to use the approximate measurement model

y≈ A(µa,,µs,)h+e, (4.16) whereµa, andµs, are our estimates for the absorption and scattering of the body. The use of incorrect realisations(µa,,µs,)in the measurement model induces modelling error [A(µa,µs)−A(µa,,µs,)]h in the model.

Since the inverse problem is ill-posed, this error may cause large artifacts to the reconstructed fluorophore image [27, 28].

4.3.1 Computation of approximation error statistics

For the computation of the error statistics, a set of samples of

x = (µa,µs,h)T(absorption, scattering and fluorophore concentrations) S= {x(ℓ), l=1, . . . , Nsamp} (4.17) were drawn fromπ(x). Two random draws ofx= (µa,µs,h)T fromπ(x) are shown in Figure 4.3. For drawing the samples of optical coefficients x, we used a Gaussian (Markov random field) prior model [15].

The samples were then used for the computation of the accurate for-ward solutions A(µ(ℓ)a ,µ′(ℓ)s )h(ℓ)and approximate forward solutions A(µa,,µs,)h(ℓ)using the diffusion approximation model (2.15)-(2.18). To obtain the samples of the approximation error ε(ℓ) = [A(µ(ℓ)a ,µ′(ℓ)s )− A(µa,,µs,)]h(ℓ), mean and covariance of εwere computed by Equations (3.20)-(3.21), Section 3.2.2.

Review of results

0 0.02 0 2 0 1

µa µs h

Figure 4.3: Two draws ofµa(first column),µs(second column) and h (third column) fromπ(x) in the 2D simulation domain.

4.3.2 Results

Figure 4.4 shows the target absorption µa (column 1) and scattering µs (column 2) parameters used in the simulations. The nominal values of absorption µa, (column 3) and scattering µs, (column 4) that were used in the estimateshcemandhaemare also shown. The true target fluorophore distributionhis shown in first column in Figure 4.5. The MAP estimates with the measurement data from the target domains in columns 1 and 2 in Figure 4.4 are shown in Figure 4.5, column 2-4, in matching order of rows. The estimates are:

(MAP-REF) The MAP-ref estimates using correct fixed (µa,µs). This corresponds to the reference estimate of conventional reconstruc-tion when (µa,µs) are known exactly (i.e, forward model used is A(µa,µs)h)).

(MAP-CEM) The MAP-CEM estimates with fixed optical coefficients are shown in column 3, Figure 4.5. This corresponds to estimate with conventional reconstruction when the nominal values of(µa,µs)are incorrect (i.e, forward model used is

A(µa,,µs,)h)).

(MAP-AEM) The MAP-AEM estimates with the same fixed optical

co-TRUE NOMINAL

µa (mm1) µs (mm1) µa,(mm1) µs,(mm1)

µa=0 µa=0.02

µ’s=0 µ’s=2

Figure 4.4: First and second columns showµa(left) andµs(right) of the bodyin the test cases 1-5 (top to bottom). Third and fourth columns show the (incorrect) absorption and scatteringµa,∗

andµs,∗that are used the in the computation of the estimates MAP-CEM and MAP-AEM.

efficients are shown in column 4, Figure 4.5. This corresponds to estimate with Bayesian approximation error model when the nomi-nal values of(µa,µs)are incorrect (i.e, forward model used is A(µa,,µs,)h).

Review of results

h(true) MAP-ref MAP-CEM MAP-AEM

0 1

Figure 4.5: First column: fluorophore distribution h of the body in the test cases 1-5 (top to bottom). The second column shows the MAP-ref estimate using the correct absorption and scattering (forward matrix A(µa,µs)). The absorption and scattering imagesµaandµsare shown in columns 1 and 2 in Figure 4.4 (rows in respective order). Third column shows the the MAP-CEM estimate using the incorrect absorption and scattering values (forward matrix A(µa,∗,µs,∗)).

Theµa,∗ andµs,∗ are shown in third and fourth column in Figure 4.4. The fourth column shows the MAP-AEM estimates using the same incorrect forward matrix A(µa,∗,µs,∗).

The relative error of the MAP estimates, Error= ∥hhtrue2

∥htrue2 ×100%, (4.18) where h is the estimated fluorophore distribution and htrue is the true fluorophore distribution are shown in Table 4.1.

Table 4.1: Relative error (%) in MAP estimates for each 2D simulation test cases (test cases are numbered from top to bottom in Figs. 1 and 2).

Case href hcem haem

1 40 40 57

2 44 64 61

3 44 133 67

4 35 116 62

5 41 154 71

As can be seen from Figs 4.4 and 4.5, the systematic errors in the conventional error model estimates (MAP-CEM) increase when the mag-nitude of the error between the true optical properties (µa,µs) and the nominal values (µa,,µs,) increase (rows from top to bottom in the fig-ures). At the same time, the approximation error approach produces estimates that are free of the artefacts and are almost as good as the refer-ence estimates (MAP-REF) which have been computed using the correct values of (µa,µs). The observation is also supported by the quantitative estimation errors in Table 4.1.

A notable feature in the reconstructions is that estimate ofhhas lower contrast in the reconstructions with the approximation error model than in the reconstructions with the conventional noise model. This can be explained by the fact that covariance of the combined noisee+εis larger than the covariance of random noise (i.e.,Γe+Γε >Γe), implying that the relative weight of the data residual term compared to the prior (or regu-larisation) term becomes smaller in the MAP estimate with the approxi-mation error model compared to the conventional noise model. Thus, the

Review of results

estimate gets, loosely speaking, drawn more strongly towards the prior mean, leading to a loss of contrast. This is also seen from the error es-timates in Table 4.1. In the first row which corresponds to the case that (µa,,µs,) are correct, the estimation error with the approximation error model is larger than with the conventional noise, and this discrepancy in the error arises from the lower contrast in the estimate of h with the approximation error model.

4.4 PUBLICATION IV: A NON-LINEAR APPROACH TO DIFFER-ENCE IMAGING IN DIFFUSE OPTICAL TOMOGRAPHY The objective in DOT difference imaging is to reconstruct change in the optical properties using measurements before and after the change. Con-ventionally, the image reconstruction is carried out using the difference of the measurements and a linearised approximation of the observation model [4, 6, 7, 13, 16, 17]. One of the main benefits of the linearised differ-ence reconstruction is that the approach has a good tolerance to modelling errors, which cancel out partially in the subtraction of the measurements.

However, a drawback of the approach is that the difference images are usually only qualitative in nature and their spatial resolution can be weak, because they rely on the global linearisation of the non-linear observation model [13].

In the fourth study, we consider a non-linear approach for difference imaging in DOT. The approach is adopted from the approach which was originally developed for EIT in [66, 67]. Using this approach, the images of the optical parameters before and after the change were reconstructed simultaneously based on the two data sets.

4.4.1 Non-linear approach for difference imaging

Consider reconstruction of δx = x2−x1, given the data y1 and y2, and models (2.10)-(2.11). Let us assume that the change in optical coefficients δx = x2x1is known to occour in ΩROI ⊆ Ω, and denote the change of optical coefficients withinΩROI byδxROI. Then, δx = MδxROI where M is an extension mappingM :ΩROI →Ωsuch that

MδxROI(r) =

δxROI(r), r ∈ROI

0 r ∈Ω\ΩROI . (4.19)

Obviously, if no ROI constraint is used, we set ΩROI = . The optical coefficients after the change x2, can now be represented as a linear com-bination of the initial state and the change as

x2 =x1+MδxROI. (4.20)

Review of results

Inserting expression (4.20) to Eq. (2.11), fixing the values of the auxiliary nuisance parametersz=z˜and concatenating the measurement vectorsy1

and y2 and the corresponding models in (2.10)-(2.11) into block vectors, leads to an observation model [66, 67]

y1

Given the model (4.21), the initial state x1 and the change δx can be si-multaneously estimated as regularisation functional of ˜x= (x1,δx)Twhich allows for separate regu-larisation functionals forx1 andδxas

f(x˜) = f1(x1) + f2(δx). (4.24) This allows the use of different models forx1andδxwhen they are known to exhibit different spatial characteristics.

4.4.2 Measurements

An experiment was carried out with the frequency domain DOT instru-ment at the Aalto University, Helsinki [56, 68]. The measureinstru-ment do-mains were cylinders with radius 35 mm and height 110 mm, see Fig.

4.6. The cylindrical phantoms corresponding to states x1 and x2 are il-lustrated in Fig. 4.6. The background optical properties were approx-imately µa,bg = 0.01mm1 and µs,bg = 1mm1 at wavelength 800 nm

for both phantoms. The cylindrical inclusions in x2, which both had a diameter of 9.5mm and a height of 9.5mm, were located such that the central planes of the inclusions coincide with the central xy-plane of the cylinder domain. The optical properties of the first inclusion were ap-proximatelyµa,inc.1=0.02mm1,µs,inc.1=1mm1 (i.e., purely absorption contrast) and the optical properties of the other inclusion were approxi-matelyµa,inc.2 =0.01mm1,µs,inc.2=2mm1(i.e., purely scatter contrast), respectively.

The source and detector configuration in the experiment consisted of 16 sources and 15 detectors arranged in interleaved order on two rings located 6 mm above and below the central xy-plane of the cylinder do-main. The locations of sources and detectors are shown with red and blue circles respectively in Fig. 4.6. The measurements were carried out at 785 nm with an optical power of 8mW and modulation frequency ω =100MHz/(2π). The log amplitude and phase shift of the transmit-ted light was recorded and the nearest measurement data from each source position were removed, leading to real valued measurement vec-torsyiR360,i=1, 2. We did not have access to measured estimate of the statistics of the noisee. In the computations we approximated ei,i= 1, 2 as zero mean with covariancesΓei, where the square root of the diagonal elements (standard deviations) of error covariances Γei were specified as 1% of the absolute values ofyi,i= 1, 2. Implying that the standard devi-ations of the measurement errors are assumed to be 1% of the measured absolute values of the log(amplitude) and phase.

4.4.3 Results

Fig. 4.7 shows 2D cross sections of 3D reconstructions obtained using the non-linear difference imaging (4.21). The ROI in this case was selected as a centre slice of height z = 22mm of the cylinder (grey patch in Fig. 4.6).

The top two rows are absorption images and the bottom two rows are scattering images. Fig. 4.8 shows the corresponding 3D reconstructions with linear difference imaging (2.14).

As it can be seen from Figs. 4.7 and 4.8 the non-linear difference

Review of results

Figure 4.6: Left: a photograph of the DOT experiment using one of the phantoms. Right: Cylin-drical phantoms x1and x2with the position of sources (red circles) and detectors (blue circles).

The grey patches on the cylinders shows the ROI.

−5

Figure 4.7: Estimates using non-linear difference imaging. The top two rows areµa and the bottom two rows areµs. First and third rows are horizontal slices of the phantom along z-axis at z = 22mm, 55mm, 99mm. The second row shows vertical slices along x-axis at x = 14mm, 35mm, 63mm. The fourth row shows vertical slices along y-axis at y = 14mm, 35mm, 63mm. The blank (white) areas are regions outside ROI.

−2

Figure 4.8: Estimates using linear difference imaging. The top two rows areµaand the bottom two rows areµs. First and third rows are horizontal slices of the phantom (height 110mm, diameter 70mm) along z-axis at z = 22mm, 55mm, 99mm. The second row shows vertical slices along x-axis at x = 14mm, 35mm, 63mm. The fourth row shows vertical slices along y-x-axis at y = 14mm, 35mm, 63mm.

Review of results

imaging shows better recovery of the inclusions when compared to the conventional linear difference imaging.

Since the conventional linear difference imaging approach is known to tolerate modelling error to some extent, the performance of this method to modelling errors was also tested in Publication IV. The results show that the presence of the modelling errors due to fixing the auxiliary nui-sance parameters z = z˜ mainly affects the estimate x1 and not δx, since x1 is the “unchanging” parameter between observationsy1 andy2 in the model (2.10)-(2.11). The non-linear difference imaging was found to tol-erate errors at least to similar extent as conventional linear difference re-constructions.

5 Summary and conclusions

In this thesis, novel computational methods for diffuse optical imaging were developed. In the first publication, the modelling of errors due to optode coupling and position errors in DOT absolute imaging was con-sidered. Errors related to optode uncertainties are known to induce severe errors in DOT due to the ill-posed nature of the problem. In this publica-tion, the Bayesian approximation error approach was used for modelling measurement system related uncertainties and it was demonstrated that it is possible to refrain from such errors using using this technique. The tolerance of the approach with different magnitude of optode errors and the sensitivity of the approach with respect to the specification of the prior model in the estimation of approximation error statistics was also tested.

In the second publication, the Bayesian approximation error approach was applied to handle errors caused by inaccurately known boundary shape. The results show that the approach could be used to estimate optical coefficients in cases where the domain shape is unavailable. The approach was also tested when both domain shape modelling and dis-cretisation error is present. The results show that the approach can handle both shape and discretisation errors simultaneously and leads to compu-tationally faster forward model.

In the third publication, modelling errors due to unknown (heteroge-neous) optical parameters in fDOT was considered. Inaccurately known optical parameters are known to cause artefacts in the reconstructed flu-orophore images in fDOT imaging when fixed incorrect values are used in the reconstructions. In this publication, the unknown optical param-eters were treated as uninteresting nuisance paramparam-eters and the mod-elling errors caused by the uncertainty in the nuisance parameters were marginalised using the Bayesian approximation error approach. The ap-proach was tested with 2D simulations and a realistic 3D mouse model.

In the fourth publication, a non-linear method for difference imag-ing that models the difference imagimag-ing problem without the linearlisa-tion approximalinearlisa-tion was proposed. The feasibility of the method was tested with simulations and experimental data from a phantom. Toler-ance to modelling errors such as domain truncation, optode coupling errors and domain shape errors was studied. It was observed that this approach estimates changes in optical coefficients with better spatial res-olution than the conventional linear difference imaging estimates. The non-linear method was found to tolerate errors at least to similar extent as conventional linear difference reconstructions.

In conclusion, novel computational modelling and reconstruction meth-ods for DOT and fDOT taking into account modelling errors were devel-oped in this thesis. The modelling errors (optode errors, unknown shape, unknown µas in fDOT etc.) considered in the thesis are some of the biggest sources of errors in diffuse optical tomography modalities that have limited their applications in clinical usage.

It was demonstrated that the application of the techniques in the pres-ence of these uncertainties improved reconstructions in DOT and fDOT imaging. We expect that these techniques can bridge the gap between laboratory setups and real life clinical situations where such model uncer-tainties are practically always present and eventually bring diffuse optical tomography closer towards clinical practises.

Few of the methods presented in the thesis have been only tested with 2D simulations. The future work flow would be to verify these methods next with realistic 3D simulations and then move on to measurements with tissue phantoms in laboratories and finally to measurements from a clinical setup. The application of the methods presented in the thesis could make absolute imaging a more practically feasible method for in-vivo studies. Also, application of the methods presented in fDOT imaging could lead to better reconstructions for in-vivo studies where the optical properties of the test subject are practically always unknown. Non-linear difference imaging has been tested with experimental phantom data in this thesis and it has shown improvement over the conventional linear

Summary and conclusions

difference imaging estimates. The next step would be to test the approach with measurements from human subjects.

Bibliography

[1] S. Arridge, “Optical Tomography in Medical Imaging,” Inv. Probl.

15, R41–R93 (1999), Topical Review.

[2] S. Arridge and J. Schotland, “Optical tomography: forward and inverse problems,”Inv. Probl.25,123010 (2009), Topical Review.

[3] A. Gibson, J. Hebden, and S. Arridge, “Recent advances in diffuse optical imaging,”Phys. Med. Biol.50,R1–R43 (2005), Topical Review.

[4] S. R. Hintz, D. A. Benaron, A. M. Siegel, A. Zourabian, D. K. Steven-son, and D. A. Boas, “Bedside functional imaging of the prema-ture infant brain during passive motor activation,” J. Perinat. Med.

29, 335–343 (2001).

[5] J. P. Culver, A. M. Siegel, J. J. Stott, and D. A. Boas, “Volumetric diffuse optical tomography of brain activity,” Opt. Lett. 28, 2061–

2063 (2003).

[6] T. N¨asi, H. M¨aki, P. Hiltunen, J. Heiskala, I. Nissil¨a, K. Kotilahti, and R. J. Ilmoniemi, “Effect of task-related extracerebral circulation on diffuse optical tomography: experimental data and simulations on the forehead,”Biomed. Opt. Express4,412–426 (2013).

[7] A. Gibson, T. Austin, N. Everdell, M. Schweiger, S. Arridge, J. Meek, J. Wyatt, D. Delpy, and J. Hebden, “Three-dimensional whole-head optical tomography of passive motor evoked responses in the neonate,” NeuroImage30,521 – 528 (2006).

[8] A. T. Eggebrecht, S. L. Ferradal, A. Robichaux-Viehoever, M. S. Has-sanpour, H. Dehghani, A. Z. Snyder, T. Hershey, and J. P. Culver,

“Mapping distributed brain function and networks with diffuse op-tical tomography,” Nat. Phot.8,448–454 (2014).

[9] B. J. Tromberg, B. W. Pogue, K. D. Paulsen, A. G. Yodh, D. A. Boas, and A. E. Cerussi, “Assessing the future of diffuse optical imag-ing technologies for breast cancer management,”Med. Phys.35,2443 (2008).

[10] D. Leff, O. Warren, L. Enfield, A. Gibson, T. Athanasiou, D. Patten, J. Hebden, G. Yang, and A. Darzi, “Diffuse optical imaging of the healthy and diseased breast: A systematic review,”Breast Cancer Res.

Tr. 108,9–22 (2008).

[11] G. Gulsen, O. Birgul, M. B. Unlu, R. Shafiiha, and O. Nalcioglu,

“Combined diffuse optical tomography (DOT) and MRI system for cancer imaging in small animals.,” Technol. Cancer Res. T. 5, 351–63 (2006).

[12] J. P. Culver, T. Durduran, D. Furuya, C. Cheung, J. H. Greenberg, and A. G. Yodh, “Diffuse optical tomography of cerebral blood flow, oxygenation and metabolism in rat during focal ischemia,” J. Cereb.

Blood Flow Metab.23, 911–924 (2003).

[13] E. M. C. Hillman, J. C. Hebden, F. E. W. Schmidt, S. R. Arridge, M. Schweiger, H. Dehghani, and D. T. Delpy, “Calibration techniques and datatype extraction for time-resolved optical tomography,” Rev.

Sci. Instrum.71,3415 (2000).

[14] M. Schweiger, I. Nissil¨a, D. A. Boas, and S. R. Arridge, “Image recon-struction in optical tomography in the presence of coupling errors,”

Appl. Opt.46,2743–2756 (2007).

[15] S. R. Arridge, J. P. Kaipio, V. Kolehmainen, M. Schweiger, E. Somer-salo, T. Tarvainen, and M. Vauhkonen, “Approximation errors and model reduction with an application in optical diffusion tomogra-phy,”Inv. Probl.22,175–195 (2006).

[16] J. C. Hebden, A. Gibson, T. Austin, R. M. Yusof, N. Everdell, D. T.

Delpy, S. R. Arridge, J. H. Meek, and J. S. Wyatt, “Imaging changes in blood volume and oxygenation in the newborn infant brain us-ing three-dimensional optical tomography,”Phys. Med. Biol.49, 1117 (2004).

[17] T. Austin, A. Gibson, G. Branco, R. Yusof, S. Arridge, J. Meek, J. Wy-att, D. Delpy, and J. Hebden, “Three dimensional optical imaging of blood volume and oxygenation in the neonatal brain,”NeuroImage 31, 1426 – 1433 (2006).

Bibliography

[18] A. Corlu, R. Choe, T. Durduran, M. A. Rosen, M. Schweiger, S. R. Ar-ridge, M. D. Schnall, and A. G. Yodh, “Three-dimensional in vivo flu-orescence diffuse optical tomography of breast cancer in humans.,”

Opt. Express15,6696–716 (2007).

[19] V. Ntziachristos, J. Ripoll, L. V. Wang, and R. Weissleder, “Looking and listening to light: the evolution of whole-body photonic imag-ing,”Nat. Biotech.23,313–320 (2005).

[20] V. Ntziachristos, C.-H. Tung, C. Bremer, and R. Weissleder, “Fluo-rescence molecular tomography resolves protease activity in vivo,”

Nat. Med.8,757–761 (2002).

[21] S. C. Davis, K. S. Samkoe, J. A. O’Hara, S. L. Gibbs-Strauss, H. L.

Payne, P. J. Hoopes, K. D. Paulsen, and B. W. Pogue, “MRI-coupled Fluorescence Tomography Quantifies {EGFR} Activity in Brain Tu-mors,”Acad. Radiol.17, 271 – 276 (2010).

[22] S. Patwardhan, S. Bloch, S. Achilefu, and J. Culver, “Time-dependent whole-body fluorescence tomography of probe bio-distributions in mice,”Opt. Express13,2564–2577 (2005).

[23] A. Koenig, L. Herv, V. Josserand, M. Berger, J. Boutet, A. Da Silva, J.-M. Dinten, P. Pelti, J.-L. Coll, and P. Rizo, “In vivo mice lung tumor follow-up with fluorescence diffuse optical tomography,” J. Biomed.

Opt. 13,011008–011008–9 (2008).

[24] S. v. d. Ven, A. Wiethoff, T. Nielsen, B. Brendel, M. v. d. Voort, R. Nachabe, M. Mark, M. Beek, L. Bakker, L. Fels, S. Elias, P. Lui-jten, and W. Mali, “A Novel Fluorescent Imaging Agent for Diffuse Optical Tomography of the Breast: First Clinical Experience in Pa-tients,”Mol. Imaging Biol.12,343–348 (2010).

[25] S. J. Erickson, S. L. Martinez, J. DeCerce, A. Romero, L. Caldera, and A. Godavarty, “Three-dimensional fluorescence tomography of hu-man breast tissues in vivo using a hand-held optical imager,” Phys.

Med. Biol.58,1563 (2013).

[26] V. Ntziachristos, A. Hielscher, A. Yodh, and B. Chance, “Diffuse op-tical tomography of highly heterogeneous media,”IEEE Trans. Med.

[26] V. Ntziachristos, A. Hielscher, A. Yodh, and B. Chance, “Diffuse op-tical tomography of highly heterogeneous media,”IEEE Trans. Med.