Mathematica in digital subtraction perfusiography of burn wounds

Iveta ((Bryjová)^a) and Stanislav ((Hledík)^b)

aDept. of Epidemiology & Public Health, Faculty of Medicine, Univ. of Ostrava, Czech Rep.
bInstitute of Physics, Faculty of Phil. & Science, Silesian University in Opava, Czech Rep.

Abstract

PIPPA (acronym for the Perfusion Image Progression and Prognosis Assessment) is a project aimed at improving the early-posttraumatic-day laser Doppler diagnostics of burn injuries. PIPPA facilitates subtraction of blood perfusion maps of a burn wound taken with an interval of a few days. Such a task of obtaining subtraction maps under clinical conditions presents a challenge since the lapse of time between the two scans makes them differ more than the algorithms for image registration can cope with. The extensive interactive and image processing capabilities of Mathematica are advantageously used in the form of Digital Subtraction Perfusiography (DSP) method.

The perfusion and/or subtraction maps overlaid over the wound area photographs permit topographic discrimination of  healing and compromised loci in burns.  The expected benefit is improving diagnostics and topographical accuracy in treatment planning, particularly in reduction of unnecessary surgery, which  is currently being clinically tested but still under development.

Skin microcirculation & LDPI

LDPI: Laser Doppler Perfusion Imaging (in vivo)

({{, , }, {Scattering on moving RBC, Δf ≡ f - (Subscript[f, 0]) = -2(Subscript[f, 0])  ((v cos θ)/c), }})

Blood perfusion diagnostic method based on Doppler frequency shift of laser beam backscattered by moving RBC

Laser beam: wavelength 630-780 nm, output power ~ 1 mW, beam diameter ~ 1 mm

Bandwidth filter 50Hz–15kHz   ⇒  maximum RBC velocity ~ 5 mm/s   ⇒  entire microvascular velocity range

For details see review by Leahy MJ et al.: Biophotonic methods in microcirculation imaging, Medical Laser Application 22 (2007),  pp 105–126

Quantitatively expressed in Perfusion Units (PU):  ∝  # of RBC passing unit area ⊥ to surface per unit time

Skin microcirculation in burns

LDPI blood perfusion maps

Medium-term (healing) period: 48 hrs up to 2–3 wks after injury: repair processes create high perfusion that can meet metabolic requirements

An important tool for estimating burn depth and reinforcing the clinical judgement of a burn wound

Kloppenberg FWH, Beerthuizen GIJM, ten Duis HJ, Burns 27 (2001), pp 359–363:
Skin microcirculation (perfusion) in burns is related to burn depth and healing time

University Hospital Ostrava: Data archives since Jan 2010 from patients with variously extensive and severe burns

ROI-averaged perfusion time evolution

Criteria for spontaneous healing within 3 wks: 
perfusion ≿ critical value  AND  (d/(dt))perfusion ≿ critical value

Box-and-whisker plot from our paper in preparation:

Call for more: maps of perfusion time evolution

Allows local assessment of perfusion temporal change, and local prognosis

Avoiding unnecessary early necrectomy locally

But! Not a straightforward problem: registration algorithms failure on vague images

The PIPPA project

Perfusion Image Progression and Prognosis Assessment

What’s the existing LDPI device output?

Four components (b–d combined in proprietary *.ldi format):

Camera photograph (separate *.jpg file)

Intensity map – shares morphological features with camera photograph

Perfusion map – most important, but much less common features

Header information (image size, resolution, etc.)

Goals

Constructing maps of local perfusion time evolution (not implemented yet by manufacturer-supplied software)

Using solely the existing LDPI modality output

Providing user-interactive environment for medical assessment

How to achieve these goals?

Choosing two scans separated typically by 48 hrs.

Embedding intensity map into photograph navigated by common features, and saving the transform.

The perfusion map transform is defined by the intensity map transform desribed above.

Registering the photographs and saving the transform.

Applying the transforms on perfusion maps.

Subtracting aligned perfusion maps.

Presentation of results in a handy and interactive way.

Not an easy task: the images to compare may  differ significantly, more than the algorithms for image registration can cope with. At this point the Mathematica’ s extensive interactive and image processing capabilities gain ground in the form of Digital Subtraction Perfusiography (DSP) method.

Implementation with Mathematica

Geometric Operations

i1=

ImageResize[i1, Scaled[0.6]]

ImageTrim[i1, {{128, 52}, {302, 265}}]

ImageRotate[i1, 15 Degree,#]&/@{Automatic,All,Full}

Manipulation and Arithmetic

i2=;

α=0.4;
ImageCompose[i1,{i2,α}]

ImageCompose[i1,{i2,α},Scaled[{.7,.3}]]

ImageDifference vs ImageSubtract

g = Image[DensityPlot[#, {x, 0, 1}, {y, 0, 1}, ColorFunction -> GrayLevel, ColorFunctionScaling -> False, PlotRange -> {0, 1}, BoundaryStyle -> None, PlotRangePadding ->None, Frame ->False, FrameTicks -> None, Mesh -> None, PerformanceGoal -> "Quality"], "Real", ColorSpace -> "Grayscale", ImageSize -> 96] & /@ {y, x}
b = Binarize[#, .5] &/@g

ImageDifference

Grayscale/Color: both image data and image are absolute value of pixel-by-pixel subtraction

gdif=ImageDifference@@g;
Row[Insert[Insert[Flatten[{g,gdif}]," - ",2]," = ",4]]
Column[Short/@(#[ImageData@gdif]&/@{First,Last})]

Binary: both image data and image are absolute value of pixel-by-pixel subtraction (logical XOR)

bdif=ImageDifference@@b;
Row[Insert[Insert[Flatten[{b,bdif}]," - ",2]," = ",4]]
Column[Short/@(#[ImageData@bdif]&/@{First,Last})]

ImageSubtract

Grayscale/Color: image data true pixel-by-pixel subtraction, but negative values clipped in image

gsub=ImageSubtract@@g;
Row[Insert[Insert[Flatten[{g,gsub}]," - ",2]," = ",4]]
Column[Short/@(#[ImageData@gsub]&/@{First,Last})]

This is exactly what we need!

Binary: both image data and image are pixel-by-pixel subtraction with -1’s clipped to 0
(logical ! ⇒)

bsub=ImageSubtract@@b;
Row[Insert[Insert[Flatten[{b,bsub}]," - ",2]," = ",4]]
Column[Short/@(#[ImageData@bsub]&/@{First,Last})]

Conclusion:

Grayscale/Color image data: ImageSubtract ≠ | ImageSubtract | = ImageDifference

Image/@{ImageData@gsub,Abs[ImageData@gsub],ImageData@gdif}

Binary image data: ImageSubtract = | ImageSubtract | ≠ ImageDifference

Image/@{ImageData@bsub,Abs[ImageData@bsub],ImageData@bdif}

Filtering

True size (64 x 64 pixel)

{#,MedianFilter[#, 2], GaussianFilter[#, {3, 0.5}],GaussianFilter[MedianFilter[#, 2], {3, 0.5}],ImageConvolve[#,SavitzkyGolayMatrix[3,2]],ImageConvolve[MedianFilter[#, 2],SavitzkyGolayMatrix[3,2]]}&[]

Factor 2 enlarged

{#,MedianFilter[#, 2], GaussianFilter[#, {3, 0.5}],GaussianFilter[MedianFilter[#, 2], {3, 0.5}],ImageConvolve[#,SavitzkyGolayMatrix[3,2]],ImageConvolve[MedianFilter[#, 2],SavitzkyGolayMatrix[3,2]]}&[]

Favorable properties of Savitzky–Golay filter:

its order (= order of the smoothing polynomial) also equals to the highest
conserved moment

Gaussian preserves 0th, 1st, but violates 2nd and higher

# of conserved moments adjustable for SG (0 = volume under, ie. total perfusion, 1 = peak mean position, 2 = peak width)

Blowing up low res LDI scan to an image

data={{0,0,0,0},{0,1,0,1},{0,0,2,0},{0,1,0,3}};
rdata=Reverse[data];
GraphicsRow[Flatten[{Style[MatrixForm[rdata],24],ArrayPlot[rdata],MatrixPlot[rdata],#[data, Mesh -> None, Frame -> None, ColorFunction -> "TemperatureMap",InterpolationOrder -> 3]&/@{ListDensityPlot,ListContourPlot}}],ImageSize->800]

Dynamic interactivity

Manipulate[ImageCompose[i1,{i2,alpha},Scaled[{tx,ty}]],
{{alpha,0,"α"},0,1},{{tx,.5},0,1},{{ty,.5},0,1}]

... and many others

Masking (morphological processing)

Interpolation

Calibration for overlays positioning

Putting it all together

The PIPPA IDE notebook with all the necessary supporting files are not available with this presentation but can only be requested at hlediks@gmail.com; expect (with apologies) the following link will not work on your computer:

Go to PIPPA notebook |->

However, you can read about the PIPPA workflow by following this link:

The PIPPA workflow |->

Conclusion

Benefits of PIPPA

Provides subtraction maps (of two perfusion maps with a time lapse between)

Allows more precise spatial and temporal analysis of blood LDI perfusion

Improvement of diagnostics & accuracy in treatment

Allows reduction of unnecessary surgery

Current status

Collaboration with manufacturer commenced

Work on Wolfram Language to C# translation: coming soon

Device development wish list (selection)

Ambient light influence elimination

Rapid scanning methods for motion artifacts elimination

Markers for scanner head position reproducibility

Scanned surface curvature correction

Software wish list (selection)

Fluid residuals in burns correction

Automatic/semiautomatic image registration

Mathematica package release (for own research purposes)

Credits

Filip Novotný (presents at this WTC) for creating the FileFormatLDPI package allowing import of native *.ldi files.

Jiří Štětinský & Hana Klosová (both Burn Center, University Hospital Ostrava, Czech Republic) for their collaboration and providing all data.

Button[Style[" Thank you for your attention. \nAny questions?",Lighter[Blue,0.6],48,LineSpacing->{1,0},TextAlignment->Center],Speak["Thank you for attention. Any questions?"],Background->Lighter[Black,1/4],Appearance->"Palette",FrameMargins->15 ]