EFFECT OF NOISE ON SIGNALS IN IN-VIVO NEAR INFRARED SPECTROSCOPY

Document Sample
EFFECT OF NOISE ON SIGNALS IN IN-VIVO NEAR INFRARED SPECTROSCOPY Powered By Docstoc
					EFFECT OF NOISE ON SIGNALS IN IN-VIVO NEAR INFRARED
                  SPECTROSCOPY




                            Alan Tan

          Supervisors: Dr St, Lawrence and Mamadou Diop

                         6 Week Project

                         April 15th 2009
                                                                Table of Contents
1.0 INTRODUCTION ...................................................................................................... 3
   1.1 HISTORY OF IN VIVO NEAR INFRARED SPECTROSCOPY .................................................................... 3
   1.2 THEORY BEHIND NEAR INFRARED SPECTROSCOPY........................................................................... 3
      1.2.2 Brain Tissue Model ..................................................................................................................... 4
2.0 METHOD .................................................................................................................... 6
   2.1 EXPERIMENTAL DESIGN ..................................................................................................................... 6
   2.2 THEORETICAL MODEL........................................................................................................................ 6
   2.3 IRF ...................................................................................................................................................... 7
   2.4 CONVOLUTION .................................................................................................................................... 8
   2.5 AMPLITUDE RATIO.............................................................................................................................. 9
   2.6 NOISE GENERATION ............................................................................................................................ 9
   2.7 LOOP ALGORITHM .............................................................................................................................. 9
3.0 RESULTS .................................................................................................................. 11
4.0 DISCUSSION ............................................................................................................ 22
   4.1 TYPES OF NOISES ................................................................................................................................23
   4.2 EFFECT OF NOISES..............................................................................................................................23
   4.3 METHOD OF NOSE REDUCTION .........................................................................................................23
   4.4 FUTURE OF IVNIRS ............................................................................................................................24
      4.4.1 Advantage of Clinical ivNIRS ...................................................................................................25
      4.4.2 Disadvantages.............................................................................................................................25
   4.5 SOURCES OF ERROR ...........................................................................................................................25
5.0 CONCLUSION ......................................................................................................... 26
6.0 COMPLETE PROGRAM CODE ........................................................................... 27
7.0 REFERENCES.......................................................................................................... 31




                                                                                                                                                                2
1.0 Introduction

        There has always been an expanding interest in the development and use of
optical devices and methods for both biological and medical application. The growth of
knowledge and understanding of fundamental optical processes has allowed the important
optical technologies to flourish with new inventions and innovations in the field.
        In vivo near infrared spectroscopy (ivNIRS) is spectroscopic method that uses the
near infrared region of the electromagnetic spectrum for both medical and pharmaceutical
diagnostics. The project dealt with ivNIRS in the non invasive assessment of brain
function through the intact human skull and its potential to be of greater use clinically in
the future.

1.1 History of In Vivo near Infrared Spectroscopy

        ivNIRS began in the field of oximetry, the study of the oxygenation of
hemoglobin in the blood. Oximetry relies on the key fact that hemoglobin can exist either
as an oxyhemoglibin or a deoxyhemoglobin, where each has its own particular absorption
characteristic, creating different peaks in spectroscopy. Initially, the work on absorption
spectra of hemoglobin derivatives was carried out in the IR and visible wavelength.
However, around the mid 1970s, Bulter and Norris interrogated a human hand using NIR
wavelengths and Jobsis reported using NIR wavelengths to record brain and myocardial
oxygen sufficiency (Rolfe, 2000). Since then, the growth, development and application of
ivNIRS has evolved rapidly.

1.2 Theory behind Near Infrared Spectroscopy.

       1.2.1 Basic Model




       Fig 1. Near-infrared beam in the basic model

        Like all spectroscopy, ivNIRS analyses the chemical content of the target sample.
In a basic model of ivNIRS, a solution containing chemicals to be analyzes is placed in a
glass curvette and positioned within the spectrophotometer (Fig 1.). Subsequently, an
interrogating near infrared electromagnetic energy is delivered to one side of the curvette.
The stream of photons from the energy propagates through the first curvette wall, through
the sample and through the second curvette wall, finally coming out with a reduced
intensity. The intensity is reduced due to the interrogating beam undergoing an
absorption by the chemical of interest. During this process, the energy of the photon


                                                                                          3
raises the molecule of the chemical to an excited state, which translates to an increase in
molecular vibration, perceived as an increase in temperature. The energy from the photon
is therefore eliminated by becoming heat. Since different chemicals respond to different
energy levels, different wavelengths elicit different level of absorption by the chemicals.
By determining how much of the intensity was lost, the following equation, which is a
derivative of the Lambert-Beer Law, can be used to determine the absorption coefficient
(Rolfe, 2000).



where,
I = reduced intensity
Io= incident intensity
ε = molar extinction coefficient (absorption coefficient)
[C] = concentration of sample
L = width of the curvette

The absorption coefficient determined pertains to the chemical only at the wavelength
and this allows oxyhemoglibin or deoxyhemoglobin to be differentiated in oximetry.

        1.2.2 Brain Tissue Model

         In clinical applications however, tissue samples instead of solutions are used most
of the time. However, this is much harder and several simplifications and assumptions
must be made and their validity has still not yet been fully evaluated. Since the sample is
not a solution enclosed in a transparent class curvette but a physically and chemically
heterogeneous tissue sample, the incident interrogating beam will have to propagate
through a much more complex structure than a simple liquid before emerging (Houten,
1996). The myriad of interfaces lead to many other effects other than simple absorption.
Particles scatter electromagnetic radiation if they have a refractive index that is different
from that of the background materials, and thus scattering happen at blood vessel walls,
collagen matrix, cell walls, and within organelles (Rolfe, 2000). As a result of this, the
distance traveled by scattered photons between the transmitter and the receiver is longer
than that traveled by the unscattered photons of the basic model (Fig 2.).




                      Fig 2. Near-infrared beam in a brain tissue model
                                                                                            4
        Therefore, it takes a much more complicated model shown below to depict the
absorption and transmission of NIR ray than equation 1.




       Where
       R = reflectance
       l = distance between transmitter and detector
       Zo = Scattering Coefficient constant (1/μs)
       t = time
       D = Scattering Coefficient constant (1/ μs)
       μa = Absorption coefficient constant

        In this investigation only the absorption coefficient is important and the other
variables have been given by the supervisors and will not be covered further
         For complex model dealing with different kinds of scatter and a longer optical
path lengths, approaches based on NIR time resolved spectroscopy (NIRTRS) have been
often used (Hoshi,2003). In NIRTRS, NIR laser pulses are input to the tissue and the
emerging intensity is detected as a function of time, providing a temporal profile of the
intensity of the transmitted light (Fig 3). Basically, it shows the population distribution of
time delay and the mean optical path length of the photons. The peak of the curve depicts
the time delay of the majority of the photons. The NIRTRS can be quantified by
estimating the absorption coefficient with fitting of the observed temporal profile to the
theoretical model of photon diffusion (Hoshi, 2003). The theoretical model is a diffusion
equation which is an approximation of the transport equation based on both the zero
boundary condition and the extended boundary condition (Diop, 2009). However, the
derivation of the theoretical is beyond the scope of the research and is not covered in this
project. An important point is that only one absorption coefficient is collected from the
ivTRS data unlike the basic model, and this absorption coefficient is the summation of all
the chromophores (HHb, O2Hb and water mainly) in the tissue and predetermined
mathematical techniques are required to extract the concentration of each chemical
component (Rolfe, 2000).

       Since different tissues and organs vary in their structure and composition, the
measurement of tissue optical properties also vary and have been an important aspect of
developing ivNIRS.




                                                                                             5
                      Fig 3. Time resolved spectroscopy measuring the frequency
                      (intensity) of photons at different time delays
2.0 Method

   2.1 Experimental Design

   1. Use Matlab to determine a program that finds a noise level acceptable during the
      time resolved spectroscopy
   2. The acceptable deviation is 10% and the maximum noise allowed yielding less
      than 5% deviation will be determined from a range of noise levels.
   3. The noise levels used will be between 5% and 50% noise at 5% interval
   4. 100 trials will be conducted per noise level, making the total number of trials
      1000 (100 x 10)
   5. Analyze the deviation with regards to the absorption coefficient
          a. Absorption coefficient is the summation of all the chromophores (HHb,
              O2Hb and water mainly) in the tissue

        Since this was a multidiscipline experiment and I was involved with the data
analysis aspect of the ivNIRS investigation, the project does not deal with the detailed
method associated with obtaining the NIRS data. The NIRS data was provided and a
mathematical program was written using MATLAB to determine the effects noise has on
the value of the absorption coefficient.
        The analysis requires both the convolution of theoretical data with the actual IRF
data as well as fitting the resulting curve to determine the absorption coefficient

       2.2 Theoretical Model

        As mentioned, the theoretical model is a solution to the photon diffusion equation
based on using the extended boundary condition. The theoretical model used pertains to
the photon diffusion of the brain. Initially, the zero boundary condition was suppose to be
used, however the program ran too slow with the condition and it was switched to the
extended boundary condition, which is a bit more complicated and is beyond the scope of
the project. The theoretical model of the brain was determined to be the following from a
primary source (Exhibit 1). The equation for the model was written on the program and a
curve was generated based on variables given by the supervisor (Image 1).

Exhibit 1: Code for theoretical model

%Theoretical Model:

v = 300; %speed of light in vaccum expressed in mm/ns
    n = 1.4; %refractive index of the medium
    c = v/n; %speed of light in the medium
    ua = 0.2; %absorption coefficient in mm-1
    us = 0.9600; %reduced scattering coefficient in mm-1
    D = 1/(3*(ua+us));%diffusion coefficient in mm
    zo = 1/(ua+us); %first point source term in mm
    zb = 2*D*(1+0.493)/(1-0.493); %2nd point source term in mm
    rho = 30; %distance from the light source on the surface in mm
    alpha = -ua*c;
    beta = c*(4*pi*D*c)^-1.5000;


                                                                                          6
     gamma = -(zo^2+rho^2)/(4*D*c);
     delta = -((zo+2*zb)^2+rho^2)/(4*D*c);
     epsilon = 1/2*(4*pi*D*c)^-1.5000*zo;
     phi = 1+2*zb/zo;

    Reflectance = 6.644e002.*(exp(alpha.*t).*(0.1180.*(beta.*t.^-
1.5000.*(exp(gamma./t)-exp(delta./t)))+0.3060.*(epsilon.*t.^-
2.5000.*(exp(gamma./t)+phi.*exp(delta./t)))));

    %Convolve the theoretical data with the IRF
ExpData = conv(IRF,Reflectance);
ExpData = ExpData(1:length(t));

       2.3 IRF

        The IRF data is a long string of numbers that represents intrinsic influences the
laser spectrometer has on the collected NIRS data. The IFR data, when graphed, will
produce a curve very similar to the theoretical model. First a subprogram was written that
pulls out in IRF file and analyzed the IRF strings and eliminated all the 0s, which should
not be included in the analysis because they represent lack of IRF activity. A command
was then written to separate the IRF data into three components so that it would take
three runs to analyze the complete IRF, and thus increasing accuracy (Exhibit 2).

Exhibit 2: Code for IFR

%Open the IRF File for Reading
[fname fpath] = uigetfile('*.asc', 'Choose Instrument Response Curve',
'C:\Documents and Settings\Mamadou Diop\My Documents\MATLAB\Data');
name = [fpath fname];
fid = fopen(name, 'rt');           % Open text file
%Read the Lines of the file
InputText=textscan(fid,'%s',10,'delimiter','\n'); % Read strings
delimited by a carriage return:'Header'
Header=InputText{1};
disp(Header);
IRF = textscan(fid,'%f',12288,'delimiter','\n'); % Read floating
numbers delimited by a carriage return:'Data'
IRF = IRF{1}; %Convert the data from cell to matrix
IRF_1 = IRF(1:4096);
IRF_2 = IRF(4097:8192);
IRF_3 = IRF(8193:12288);
IRF = (IRF_1 + IRF_2 + IRF_3)/3;

timeRes = 3.010e-003; %Time/channel

%Finds start and end points of IRF and creates truncated time
counterVal = 0;
counter = 1;
while counterVal == 0
    counterVal = IRF(counter);
    counter = counter + 1;
end




                                                                                         7
endpoint = counter + 10/timeRes; %Only take the 10 ns from the
bigining
IRF = IRF((counter + 1):endpoint);

%Time Series
%Time Series for data and IRF
t = counter + 1:endpoint;
t = t';
t = t*timeRes; %Convert the time into ns

       2.4 Convolution

       This step involved convoluting the theoretical model with the IRF data. Finally, a
convolution function was used to convolute the theoretical brain diffusion model with the
IRF data. By convoluting the theoretical model with the IRF data, an actual experimental
model of what the data should be if the brain was analyzed by the ivNIRS was generated,
but without the noise.

Exhibit 3: Code for Convolution

%Convolve the theoretical data with the IRF
ExpData = conv(IRF,Reflectance);
ExpData = ExpData(1:length(t));
plot(t, ExpData, 'r', t, IRF, 'b')




Image 1. Comparison between the Theoretical and the IRF curve (Red = IRF, Blue = Theoretical
Data)




                                                                                           8
       2.5 Amplitude Ratio

        The IRF data and the theoretical model produce a similar curve, but their sizes are
very dissimilar where the curve of the IRF has higher amplitude than that of the
theoretical model. As a result, an amplitude multiplier was determined that would make
the time resolved spectrometry curve of the theoretical model as close to that of the IRF
data as possible. In this case, it was determined to be 10^22 by trial and error.

       2.6 Noise Generation

       This stepped involved writing a random noise generator subprogram that would
generate random noises at 10 different noise levels (Exhibit 4). The random noise
generator generated a Poisson random number for every 0.001 increment in time for each
convoluted curve. For each noise level, only a percentage of the each random number
corresponding to the noise level was taken, for example, if the number 4 was generated at
5% noise level, the actual number taken will be 5% of 4 or 0.02. In addition, each
generated random noise had to be divided by the amplitude determined in the step above
to compensate for the small amplitude of the theoretical model. When the amplitude was
not compensated for, the noise generated was too large and masking it on top of the
curves created distortions that completely drowned the original curve. The noise
generated was then masked on top of the convoluted curve to make an actual
experimental model that includes the noise

Exhibit 4: Code for noise Generation for one loop
%Noise
noise = 0.05./(1e22).*random('poisson',1,size(t));
%5% is the noise level here and 1e22 is the amplitude ratio (10^22)
%noise = 0.20*randn(size(t));
%'Normalize' IRF
Data = ExpData + noise;

       2.7 Loop Algorithm

        The above steps had to be repeated at 10 different noise levels at 100 trials per
level. As a result, a loop algorithm was written that ran the program 10 times
incrementally between 5% and 50% noise level (Exhibit 5). Furthermore, within the noise
loop, a trial loop was written that ran the convolution, noise generation and noise
masking 100 times for each noise level. Therefore, 1000 (100 x 10) curves were
generated following all the loops. Image 2 shows the effect of 10 noise levels on the
convoluted curve. The distortion increases with the noise level.

Exhibit 5: Loop Algorithm for curve fitting for 100 loops

%Fitting

     %Start and Endpoints of fit:

     [PeakValue PeakPosition] = max(ExpData);



                                                                                          9
    spt = PeakPosition;
    sptVal = PeakValue;
    while sptVal > 0.01*PeakValue
        sptVal = ExpData(spt);
        spt = spt - 1;
    end


    %Find endpoint for fminsearch
    ept = PeakPosition;
    eptVal = PeakValue;
    while eptVal > 0.01*PeakValue
        eptVal = ExpData(ept);
        ept = ept + 1;
    end
%100 repetitions for one noise level counted by ct
for ct=1:100
    %Noise
    noise = 0.05./(1e22).*random('poisson',1,size(t));
    %noise = 0.20*randn(size(t));

Data = ExpData + noise;

%Fitting the Data and getting an estimate for the absorption
coefficient
    Starting = [200 0.08 1];

options=optimset('Display','final','MaxFunEvals',100000,'TolFun',1.0000
e-07,'TolX',1.0000e-07,'MaxIter',100000);
    Estimates = fminsearch(@TPSF2param_fit, Starting, options, t, Data,
IRF, spt, ept);
    histua(ct) = Estimates(2);
end




      Curve Fitting and Distribution curve


     Image 2: Effect of 5% to 50% noise (from top left to bottom right respectively) on
     the convoluted curve                                                                 10
        Once the experimental model was established with the noise 1000 times, each
curve was fitted to determine the estimate of the absorption coefficient (Exhibit 5). The
program was written to fit the theoretical model onto the experimental model and
calculate the absorption coefficient. The absorption coefficient determined for each noise
level was then programmed to be displayed as a histogram using 100 increments of
distribution (bin = 100). Finally, the histogram was converted using curve function to a
curve and statistical tests were carried out using mat lab functions.


3.0 Results




Graph 1. Distribution of absorption coefficient with 5% noise




                                                                                        11
Graph 2. Distribution of absorption coefficient with 10% noise




Graph 3. Distribution of absorption coefficient with 15% noise




                                                                 12
Graph 4. Distribution of absorption coefficient with 20% noise




Graph 5. Distribution of absorption coefficient with 25% noise




                                                                 13
Graph 6. Distribution of absorption coefficient with 30% noise




Graph 7. Distribution of absorption coefficient with 35% noise




                                                                 14
Graph 8. Distribution of absorption coefficient with 40% noise




Graph 9. Distribution of absorption coefficient with 45% noise




                                                                 15
Graph 10. Distribution of absorption coefficient with 50% noise

         Graphs 1 to 10 show the histogram of the value of absorption coefficients
collected between 5% to 50% noise levels. The distribution of the histogram indicates
that most of the absorption coefficients collected fall within the average absorption
coefficient value of 0.2μm-1. As the value of absorption coefficient deviates from this
average, the frequency or the intensity of the absorption coefficients decreases with
greater deviation. However, as shown by the gradual change from 5% to 50% noise
level, the deviation from the average absorption coefficient increases. A more in depth
analysis of the trend will be done with the normal curves constructed from the data
collected in the histogram.




                                                                                          16
Table 1: Standard deviation and mean of absorption coefficient collected between 5%
and 50% noise level
                                          Standard Deviation (μm-
Noise Level          Mean (μm-1)          1)                         % Deviation
                5%                  0.2                    0.0137                    6.85
               10%                  0.2                    0.0567                   28.35
               15%                  0.2                    0.1791                   89.55
               20%                  0.2                    0.2699                  134.95
               25%                  0.2                     0.311                   155.5
               30%                  0.2                    0.3958                   197.9
               35%                  0.2                    0.4024                   201.2
               40%                  0.2                    0.4396                   219.8
               45%                  0.2                    0.4768                   238.4
               50%                  0.2                      0.491                  245.5

Table 1 depicts the gradual increase in standard deviation with increasing noise while the
mean remains constant for all of the noise levels. From this we determine that for the
deviation from the average absorption coefficient to be less than 10%, the noise level has
of the spectrometer has to be less than 10%. Even at 10% noise, there is 28% deviation,
so the maximum tolerated noise level lies somewhere between 5% and 10% noise level.




Graph 11. Distribution of absorption coefficient with 5% noise




                                                                                            17
Graph 12. Distribution of absorption coefficient with 10% noise




Graph 13. Distribution of absorption coefficient with 15% noise




                                                                  18
Graph 14. Distribution of absorption coefficient with 20% noise




Graph 15. Distribution of absorption coefficient with 25% noise




                                                                  19
Graph 16. Distribution of absorption coefficient with 30% noise




Graph 17. Distribution of absorption coefficient with 35% noise


                                                                  20
Graph 18. Distribution of absorption coefficient with 40% noise




Graph 19. Distribution of absorption coefficient with 45% noise


                                                                  21
Graph 20. Distribution of absorption coefficient with 50% noise

         The histograms by themselves were hard to analyze. Therefore, they were
converted into a normal distribution curve as shown from graph 11 to 20. However the
normal curve is not a direct conversion from the histogram as that is beyond the scope of
what is being investigated. The normal distribution is merely a representation of the
approximate distribution of the absorption coefficient at each noise level derived from the
mean and standard deviation of the histogram at each noise level. As seen from the
distribution, although the curve by itself does not show it, the width of the curve
increases with the noise level, demonstrating the gradual increase in standard deviation.
For a 5% noise level, the range of deviation is between 0.1 and 0.3. However, for a 50%
noise level, the range of deviation has increased to between-2 and 2. In addition, with the
increase of noise level, the frequency, or intensity of the absorption coefficient close to
the mean decreases as more and more absorption coefficient values deviate from the
mean. For a 5% noise level, the frequency of the absorption coefficients that are 0.2μm-1
is 30. However, for a 50% noise level, the frequency of the absorption coefficients with
this mean value has decreased to 0.85 (Please note that this is all relative and is not the
actual frequency depicted in the histograms)

4.0 Discussion

        The objective of this project was to determine the effect of noise on the absorption
coefficient of a brain time resolved ivNIRS. An increase in noise level increased the



                                                                                         22
fluctuation of the time resolved curve and contributed to more deviation in the estimated
absorption coefficient value.
        The noise of the laser is inherent and cannot completely be avoided in ivNIRS due
to various influences of quantum noise and fluctuation of various technical origins. The
quantum noise is intrinsic to all forms of light associated with the spontaneous emission
in the gain medium. The technical noise arises from excess noise of the pump source,
from vibrations of the laser resonator, or from temperature fluctuations and can be
controlled based on the specification of the laser spectrometer (Coldren,1995).

      4.1 Types of noises
      There are multiple kinds of laser noises in ivNIRS that amalgamate to create a
combined noise effect.

        Intensity noise
        Intensity noise is strictly the noise of the optical power. Intensity noise of a laser
results from the quantum noise, associated wit laser gain and resonator losses, and partly
from technical noise sources such as excess noise of the pump source, vibrations or
resonator mirrors, and thermal fluctuations in the gain medium (Yarviv,1997). As
mentioned, this intrinsic noise cannot be avoided but it can be reduced however by
changing the operation condition, as the noise becomes weaker at high pump powers,
where relaxation oscillations are strongly damped.

       Phase noise
       The output of a single frequency laser is not exactly monochromatic but exhibits
some sort of phase noise (Yarviv,1997). The origin of this noise is caused by spontaneous
emission of the gain medium, and quantum noise associated with optical losses. In
addition, there can be technical noise influences, e.g. due to vibrations of the cavity
mirrors or to temperature fluctuations (Yarviv,1997). The noise is responsible for
continuous frequency drift, or as sudden phase jumps, or a combination of both.

        Frequency noise
        The term frequency noise refers to random fluctuations of the frequency of an
oscillating signal and is directly related to the phase noise (Yarviv,1997). Phase noise and
frequency noise are established as two different ways describing the same phenomenon.

        4.2 Effect of noises
        The most obvious effect of laser noise in this project is the increase in deviation in
the absorption coefficient due to noise. However, laser noise also limits the data
transmission rate achievable with the laser spectrometer, leading to slower data
collection. In addition, it increases beam pointing fluctuations and pulse energy variation
fluctuation which leads to the low precision of the NIRS data (Coldren,1995).

        4.3 Method of Nose Reduction
        The both the quantum and technical noise can be reduced to an extent to enhance
the accuracy of the NIRS data. In general, quantum noise can be reduced by increasing
the intracavity power level and by minimizing optical losses. Therefore, theoretically, the



                                                                                            23
signal to noise ratio can keep increasing with greater laser power level. However,
clinically, the NIR laser is exposed to live subjects and a safety standard is in place that
limits the intensity of the laser you can input in the subject as dictated by the NSERC
(Diop, 2009). In addition, technical noise influences can be reduced by temperature
stabilization of the setup, or by using a low-noise pump source. Furthermore, laser
parameters can be optimized so that the laser reacts less strongly to certain noise
influences. Nevertheless, theses tweaks and adjustments require upgrades to the hardware
and software that can prove to be very costly.
        Apart from the aforementioned slight adjustments, an easier method of noise
reduction is laser stabilization. Laser stabilization schemes usually involve some kind of
electronic feedback system, where fluctuations of some parameters are converted to an
electronic signal, which is feedback to act on the laser (Coldren,1995).




       Figure 1: Diode-pumped solid-state laser with a feedback system stabilizing the
output power.

        The output power of a laser may be stabilized with a scheme as shown in
Figure 4. The laser power is monitored with a photodiode and corrected e.g. via control
of the pump power or the losses in or outside the laser resonator . In this way, the
intensity noise under steady-state conditions can be reduced. It is also possible to reduce
intensity noise by acting on the output beam instead of the laser itself and is called noise
eating (Coldren,1995). The stability which is achieved with the active systems is
determined by factors such photodetection noise, the bandwidth of control elements, the
design of the feedback electronics, and the stability of the reference standards
(Coldren,1995).

       4.4 Future of ivNIRS

        The importance of NIR time resolved spectroscopy of the brain lies in its potential
to be of wide clinical use in the near future. As mentioned, ivNIRS is very applicable in
oximetry, and is capable of analyzing the oxygen content of hemoglobin in the brain and
determine the HHb, O2Hb and oxidized cytochrome aa3 concentration of the test area.
By employing several wavelengths and time resolved methods, blood flow, volume and
oxygenation can be quantified. Therefore, in clinical scenarios, changes in blood
hemoglobin concentrations from the norm can be determined through ivNIRS.
Applications of oximetry by NIRS methods include the detection of illnesses which affect
the blood circulation, the detection and assessment of breast tumors, and the optimization
of training in sports medicine.


                                                                                           24
       4.4.1 Advantage of Clinical ivNIRS

        There are many advantages to implementing ivNIRS in clinical studies. A major
advantage is its very economical. The cost of a simple NIR spectrometer is only 10000
dollars and is perfect for quick oximetry diagnosis of the brain to see if the scatter and
absorption is within the normal range. Other devices such as MRIs are substantially more
expensive to operate and usually produce a long waiting list of patient who may have to
wait for weeks before getting a scan. Furthermore, ivNIRS requires little or no time to
retrieve the results. The NIR light travels from the transmitter to the receiver in
picoseconds and repetitions can be completed within a very short time. The NIR machine
is also very portable and can have a dimension of 10x10x10, allowing the equipment to
go to the patient while other imagining equipments are extremely bulky and requires the
patient to go to the equipment. In summary, ivNIRS’ cost-efficiency, portability and fast
scan time makes it ideal for quick preliminary diagnoses in clinical settings. Therefore,
ivNIRS will prove very useful in carrying out fast and simple diagnosis of brain
injury/disease in the ICU or ER where efficiency and flexibility is paramount in saving
the lives of the patients.

       4.4.2 Disadvantages

        Although ivNIRS has many advantages, it has its share of disadvantages. One of
its biggest disadvantages is its relatively weak penetration. With a conventional ivNIRS
under the NSERC standard, the NIR radiation can only penetrate 3cm into the skull.
Since child skulls are relatively thick, ivNIRS does not pose a problem analyzing the
brain of children. However, this poses a problem in analyzing adult brains since the
radiation might not even make it to the cerebral cortex due to the thickness of the adult
skull, and a clear data of the brain cannot be collected. In addition, the data, whether a
curve or spectrum, collected from a ivNIRS cannot be interpreted immediately since the
data and the absorption coefficient determined is a summation of all the chemicals
involved. Different predetermined mathematical models are required to fit the curve
before any information about oxygenation, chemical composition, or abnormality can be
quantified and qualified.

         The photons of the ivNIRS are assumed in the model to exhibit a mix of three
scattering patterns at the brain as shown in Fig 2. The first pattern is simply the photon
passing through the tissue unscattered. The second pattern involved the photon scattering
but still maintaining the same net direction as the input beam. The final pattern depicts
the photon scattering to the point where they emerge with a net direction of 180 degrees
to the.

       4.5 Sources of Error

        Most of the error comes from the ability to use the program Matlab. Due to the
complex functions and difficult codes that need to be learned, the investigation was
simplified and generalized.



                                                                                         25
        First of all, the random noise generated is a “pseudo random”, meaning that the
random number is not “truly random”., but is dictated by an algorithm. Therefore, every
time you run the program, the random noises generated are the same random noises every
single time per run. Therefore, the random noises among each of the 1000 trials will be
different, but if the program was run again, each of the 1000 trials will have the same
random noise assigned to it as the first run.
        In addition, the program is not optimized. Since there are many codes that are not
written optimally, the run speed of the program is very slow. It takes 10 hours to run the
program once due to the convolution and the fitting being done 1000 times. For the
program used in this investigation, the average convolution/fitting time for one trial is
~45 seconds. This is too long and in the future, better knowledge of Matlab will enable
the creation of a better program with faster run time using more advanced codes and
functions that are less time consuming than the fundamental ones.
        Last but not least, as mentioned before, the normal curve generated is not a fitted
curve to the histogram, and is not a real representation of the actual data. Better
understanding of and prowess with Matlab I the future will make these conversions
possible and yield actual data instead of representative ones.

       5.0 Conclusion

         ivNIRS is a discipline of that has a huge potential in the future. In this
investigation, the noise associated with the spectrometer was studied and it was
determined that maintaining a low noise and a high SNR is crucial to the accuracy of the
data obtained. Although much research into ivNIRS is still required for it to be a viable
tool in the medical and clinical endeavor, the demand for newer technologies and better
understanding is forever insatiable and the demand for more convenient, safe, and
economical clinical tools to improve the care of patients is never ending.




                                                                                        26
6.0 Complete Program Code

clear all;
close all;

histua = zeros(1,100);

%Program to fit Analytical Model TPSF for Time-Domain Reflectance
convolved
%with a real IRF

%Open the IRF File for Reading
[fname fpath] = uigetfile('*.asc', 'Choose Instrument Response Curve',
'C:\Documents and Settings\Mamadou Diop\My Documents\MATLAB\Data');
name = [fpath fname];
fid = fopen(name, 'rt');           % Open text file
%Read the Lines of the file
InputText=textscan(fid,'%s',10,'delimiter','\n'); % Read strings
delimited by a carriage return:'Header'
Header=InputText{1};
disp(Header);
IRF = textscan(fid,'%f',12288,'delimiter','\n'); % Read floating
numbers delimited by a carriage return:'Data'
IRF = IRF{1}; %Convert the data from cell to matrix
IRF_1 = IRF(1:4096);
IRF_2 = IRF(4097:8192);
IRF_3 = IRF(8193:12288);
IRF = (IRF_1 + IRF_2 + IRF_3)/3;

timeRes = 3.010e-003; %Time/channel

%Finds start and end points of IRF and creates truncated time
counterVal = 0;
counter = 1;
while counterVal == 0
    counterVal = IRF(counter);
    counter = counter + 1;
end

endpoint = counter + 10/timeRes; %Only take the 10 ns from the
bigining
IRF = IRF((counter + 1):endpoint);

%Time Series
%Time Series for data and IRF
t = counter + 1:endpoint;
t = t';
t = t*timeRes; %Convert the time into ns

%Theoretical Model:

v = 300; %speed of light in vaccum expressed in mm/ns
    n = 1.4; %refractive index of the medium
    c = v/n; %speed of light in the medium



                                                                         27
    ua = 0.2; %absorption coefficient in mm-1
    us = 0.9600; %reduced scattering coefficient in mm-1
    D = 1/(3*(ua+us));%diffusion coefficient in mm
    zo = 1/(ua+us); %first point source term in mm
    zb = 2*D*(1+0.493)/(1-0.493); %2nd point source term in mm
    rho = 30; %distance from the light source on the surface in mm
    alpha = -ua*c;
    beta = c*(4*pi*D*c)^-1.5000;
    gamma = -(zo^2+rho^2)/(4*D*c);
    delta = -((zo+2*zb)^2+rho^2)/(4*D*c);
    epsilon = 1/2*(4*pi*D*c)^-1.5000*zo;
    phi = 1+2*zb/zo;

    Reflectance = 6.644e25.*(exp(alpha.*t).*(0.1180.*(beta.*t.^-
1.5000.*(exp(gamma./t)-exp(delta./t)))+0.3060.*(epsilon.*t.^-
2.5000.*(exp(gamma./t)+phi.*exp(delta./t)))));

    %Convolve the theoretical data with the IRF
ExpData = conv(IRF,Reflectance);
ExpData = ExpData(1:length(t));

plot(t, ExpData, 'r', t, IRF, 'b')
xlabel('Time (ps)')
ylabel('Intensity')

%Fitting

    %Start and Endpoints of fit:

    [PeakValue PeakPosition] = max(ExpData);
    spt = PeakPosition;
    sptVal = PeakValue;
    while sptVal > 0.01*PeakValue
        sptVal = ExpData(spt);
        spt = spt - 1;
    end


    %Find endpoint for fminsearch
    ept = PeakPosition;
    eptVal = PeakValue;
    while eptVal > 0.01*PeakValue
        eptVal = ExpData(ept);
        ept = ept + 1;
    end

for noisefraction = 0.05:0.05:0.5;
for ct=1:100
    %Noise
    noise = noisefraction./(1e22).*random('poisson',1,size(t));
    %noise = 0.20*randn(size(t));

%'Normalize' IRF
%IRF = IRF - mean(IRF(1:100));
%Data = Data - mean(Data(1:100));


                                                                     28
Data = ExpData + noise;

%Fit

       Starting = [200 0.08 1];

options=optimset('Display','final','MaxFunEvals',100000,'TolFun',1.0000
e-07,'TolX',1.0000e-07,'MaxIter',100000);
    Estimates = fminsearch(@TPSF2param_fit, Starting, options, t, Data,
IRF, spt, ept);
    histua(ct) = Estimates(2);
end
end
hist(histua);

%end

%Fit Function

%IRF convolved with a theoretical Reflectance function
function sse=biGvFit(params,Input,Actual_Output,IRF,spt,ept)
amp = params(1);
ua = params(2);
us = params(3);

       v = 300; %speed of light in vaccum expressed in mm/ns
       n = 1.4; %refractive index of the medium
       c = v/n; %speed of light in the medium
       D = 1/(3*(ua+us));%diffusion coefficient in mm
       zo = 1/(ua+us); %first point source term in mm
       zb = 2*D*(1+0.493)/(1-0.493); %2nd point source term in mm
       rho = 30; %distance from the light source on the surface in mm
       alpha = -ua*c;
       beta = c*(4*pi*D*c)^-1.5000;
       gamma = -(zo^2+rho^2)/(4*D*c);
       delta = -((zo+2*zb)^2+rho^2)/(4*D*c);
       epsilon = 1/2*(4*pi*D*c)^-1.5000*zo;
       phi = 1+2*zb/zo;


Fitted_Curve = conv(IRF,
amp.*exp(alpha.*Input).*(0.118.*(beta.*Input.^-
1.5000.*(exp(gamma./Input)-
exp(delta./Input)))+0.306.*(epsilon.*Input.^-
2.5000.*(exp(gamma./Input)+phi.*exp(delta./Input)))));
%Fitted_Curve = interp1(t_temp, Fitted_Curve_pre, Input);
Fitted_Curve = Fitted_Curve(1:length(Input));
Error_Vector = Fitted_Curve(spt:end) - Actual_Output(spt:end);
sse=sum(Error_Vector.^2);

% Graphing Function

pdfMean = 0.2;
pdfStd = 0.4910;
pdfNormalSS = randn(1,100) * pdfStd + pdfMean;


                                                                        29
xRange = floor(min(pdfNormalSS)):0.1:ceil(max(pdfNormalSS));
pdfNormalFitted = normpdf(xRange, pdfMean, pdfStd);
 plot(xRange, pdfNormalFitted);
xlabel('ua (x10^{-3}m^{-1})')
ylabel('Intensity')




                                                               30
7.0 References

Rolfe, P. (2000). In Vivo Near-Infrared Spectroscopy. Journal of Biomedical
       Engineering, 2, 715-754.

Hoshi, Y. (2003). Functional near-infrared optical imaging:Utility and limitations in
       human brain mapping . Journal of Physical Physiology, 40, 511-520.

Coldren, L. A. (1995). Diode Lasers and Photonic Integrated Circuits . New York: John
       Wiley.

Yarviv, A. (1997). Optical Electronics in Modern Communications. New York: Oxford
       University Press.

Houten, V., P, J., Bernado, D., & Stevenson, D. (1996). Imaging Brain Injury Using
      Time-Resolved Near Infrared Light Scanning. . Pediatric research, 39(3), 470-
      476. Retrieved March 2, 2009, from
      http://www.pedresearch.org/pt/re/pedresearch/abstract.00006450-199603000-




                                                                                        31

				
DOCUMENT INFO