Appunti Medical Support System for Chronic Diseases - Prof. Lanzarone

MEDICAL SUPPORT

SYSTEM FOR

CHRONIC DISEASES

Prof. Ing. Ettore Lanzarone

Master Degree in Engineering and Management for Health

Academic Year 2018/2019

Index

1. Introduction

a. Estimation, Stochastic approach

b. Dynamic system, differential approach

c. Decision, optimization approach

d. Validation process

2. Lumped Parameter Model

3. Haemodialysis

4. Bioimaging

a. Diffusion-Weighted MRI, the IVIM Model

b. Computed Tomography Angiography, aortic stiffness

c. Mock loop circuits for in vitro studies of cardiovascular system

5. Artificial Pancreas

6. ECMO

7. Monte Carlo 27/02/2019

Introduction

1. Hemodialysis

There are some pathologies caused by haemodialysis, in particular due to the therapy discontinuity,

so the idea is to develop some mathematical tools so that we can reduce intra-dialysis

complications, and only then to customize the tool for each patient.

2. Lumped parameter models of the cardiovascular system

Fluid-dynamic behaviour of blood can be characterized in terms of a set of lumped parameters in

each segment; the resulting mathematical description is equivalent to that of an electric circuit, in

this way we can easily study:

a. the impact of alterations in model parameters

b. the interaction of the system with external devices

3. Bioimaging

We have a lot of images taken at different times, so the idea is to develop a tool that’s a marker for

some pathologies and can say something about potential pathologies

4. Aortic compliance evaluation

The idea is estimate functional parameters from a set of images acquired over the cardiac cycle on

a patient, so:

a. identification the radius profile in several sections

b. mechanical model to get the compliance given the radius profile

5. Extra Corporeal Membrane Oxygenation

The goal is to provide the best set of parameters to the ECMO machine and alarms if the conditions

become critical.

6. Artificial Pancreas

The goal is to perform an automatic blood glucose control, it requires an estimation of the patient’s

need and his conditions, together with a decision model to set the released quantity based on the

estimated needs Estimation, stochastic approaches

Let us assume to have a model characterized by a set of parameters, how to estimate the model

parameter from data? How to validate and compare different parameters value? How to make

statements about estimated values and assumptions?

Moment estimation

With this approach we use the observation to determine the moments of the underlying parameter

distribution (which is assumed to be known).

What we have: known statistical distribution and observations

What we don’t have: parameters of this distribution

How to estimate them: with empirical moments $

The model assumes that observation X , X , … come from a normal distribution N(, ), the

1 2

parameters of the density are estimated by fitting the empirical moments, so empirical moments

determine the parameters of the distribution, in this example we have a gaussian distribution

defined by to parameters, so we need two empirical moments that can be mean and variance.

Regression

We assume a function structure that links two variables X and Y for which we have a set of

observations (X , Y ), (X ,Y ), and so on (X and Y can be either single variables or vectors of variables).

1 1 2 2

Under the assumption for the function we aim at estimating the coefficients of the

: = ()

function. What we have or assume: the general function

What we don’t have: parameters of the function

How to estimate: minimization of error

2

How to evaluate the goodness: R e AIC

For each observation we have: , the best coefficient a is the one that minimizes

= + +

, , ,

the overall error that is the sum of the squared errors: $

{Σ }

( − − )

, , ,

Σ ( − ̅ )

, , ,

: = $

Σ ( − ̅ )

9 , ,

> :̅

= ? −

The same idea can be applied both to multiple linear regression (with several input variables) and

to non-linear functions. $

To evaluate the goodness of coefficients we use the parameter typically values greater

[0,1],

than 0.9 denote a good matching.

Stochastic approaches ?) $

Σ ( −

, ,

$

= ?)

,∗ $

Σ ( −

,

2

Why can R be smaller than 1? NOT RIGHT FUNCTION ONE STRANGE POINT

NOISY OBSERVATIONS ASSOCIATED WITH A

NO – CHANGE MODEL HIGH ERROR

OK – KEEP MODEL OK OR NO DEPENDING

ON THE CASE

The Q-Q plot is the graphical representation of the

quantiles of a specific distribution. It compares the

cumulated distribution of the observed variable with

the normal cumulated distribution

In case of multiple linear regression with several input variables to evaluate the goodness of the

model we can use the Aikake Information Criterion, it associates a score for a given combination of

included input variables and it chooses the combination with the lowest score.

Likelihood function

We consider independent and identically distributed (iid) observations coming from the same

density function with unknown parameters.

What we have: assumed density function and iid observations

What we don’t have: parameters

How to estimate them: calculating the likelihood

The likelihood is defined by the product of the density function over the iid observation; in each

factor of the product the value is known, thus the likelihood function L depends on the unknown

,

model parameters, in the Poisson case the unknown parameters is so we have only one

,

parameter: MN O MN O MN O

P S T

( , , … , ; ) = × × …×

H $ J ! ! !

H $ J

MN O

() ∝ W X

, ∑ O

X X

J JO?

∑

O O O O

∙ ∙ … ∙ = = =

J

P S T X X

JO?

MJN

() ∝

Sometimes, instead of the likelihood function its logarithm is considered (if we apply the

,

logarithm operator to a monotone function, we don’t change the maximum and minimum position).

For the Poisson case we have: JO?

MJN

() = ?

ℓ() = − + log

The best parameters that we can estimate are those that maximize the likelihood function or

,

equivalently, the logarithm likelihood function ℓ ?

ℓ = − +

?

̂ =

For the Normal case we have: S

1 (dMN)

M

() = S

$e

√2 S S S

1 1 1

(d (d (d

MN) MN) MN)

T

P S

M M M

( , , … , ; ; ) = ∙ ∙ …∙

S S S

$e $e $e

H $ J √2 √2 √2

S S

1 (d (d

MN) MN) H

X X S

∑ (d

M M M MN)

MH MJ X X

(; ) = W ∝W =

S S S

$e $e $e

√2

, ,

(? 2?

$ $

− + )

ℓ(; ) = − − $

2

For a function minimization with two variables we can use the Hessian matrix that automatically

gives a maximum in the point with null gradient, this is called Maximum Likelihood Estimator (MLE).

The MLE can be applied when observations depend on other observations, for example in the case

of a dynamic system.

Example: $ )

~( + ,

The likelihood is the product of the conditional densities evaluated in each observation, so each

term of the product is: S S

1 1

(k MN ) (k MlO Mm)

X X X X

=

S S

e e

√2 √2

The maximization is performed with respect to parameters a and b.

In this case are not iid but they are conditionally to , so it’s possible to compute the likelihood

, ,

as their product. Once again we have to make assumptions on the shape of density and relationship

between X and Y

This approach can be used to estimate the coefficients of differential models based on the

n (,

observation of state variables. Let’s consider a system characterized by: where f is a

= ),

o

given function whose set of coefficients is under the assumption that we have observation at

,

equally spaced discrete time intervals the system can be discretized with time step

∆, ∆:

−

,sH , ( )

= ,

o , ,

∆ ( )∆

= + ,

,sH , o , ,

Then, we include an error associated to each observation:

( )∆

= + , +

,sH , o , , ,sH

A typical choice is to consider a Gaussian White Noise for errors (each error follows an independent

$

Gaussian density with null mean value and variance ), when adopting such noise, the resulting

stochastic process is named Wiener process. Thus, the conditional densities are:

$

( )∆,

~( + , )

,sH , o , ,

The likelihood is the product of the conditional densities at each observation, which is function of

the set of coefficients once again the maximization is performed with respect to these

,

coefficients.

Bayesian estimation What we have: prior density and density function

What we don’t have: parameters

How to estimate: with the Bayesian theorem

It’s based on the Bayes’ theorem that gives a method to modify the trust level of a given hypothesis,

in light of a new information. Denoting with the null hypothesis, that can be assumed or deduced

u

by previous observations, and with the empirical data observed:

)(

(| )

u u

|)

( =

u ()

This is the conditioned density, that is the probability that occurs once we know has occurred

u

is the prior probability, it’s the probability distribution of the parameter of interest

( )

u

representing the probability to have certain possible values of this parameter without considering

the sampled observed data.

) is the likelihood function, it’s the probability distribution we can give to sampled observed

(|

u

data in correspondence with each specific value that the parameter of interest can assume. It

represents the probability to have some observed data, knowing the parameter of interest.

|) is the posterior probability, it’s the probability distribution of the parameter of interest

(

u

and it represent the probability to have some possible value of this parameter after having looked

at the sampled observed data.

If both the prior and the likelihood distribution support a certain value of , then we can be sure

u

that also the posterior distribution support it. But if a certain value of is not supported by the

u

Stochastic approaches

prior or by the likelihood function or by both, we can be sure that also the posterior distribution

doesn’t support it. We must remember that:

)( )

(| = (| ∩ )

u u u u

.It’s applied considering that and represent the set of parameters to be estimates and the

It is applied considering that A and B represent the data and the model

u

observed data and the structure of the model is assumed (density to write the likelihood),

parameters to be estimated:

removing the denominator we have: (|) ∝ ()(|)

With respect to MLE, we have some advantages and disadvantages:

• denotes the set of parameters to be estimated

• D denotes the data (set of observations)

• it allows to determine the entire probability density function of model parameters

• The structure of the model is assumed (density to write the likelihood)

• it allows to formalize and integrate the previous knowledge of parameters

• it’s sensible to prior density, so prior knowledge has to be proper managed

• it can be applied when each observation depends on the value at the previous time instant,

so when MLE can be applied also the Bayesian approach can be applied and vice versa.

Dynamic systems, differential approaches

(J)

n nn

(), (), ()

An Ordinary Differential Equation (ODE) is the goal is to

{, (), … , | = 0;

find the solution y(t) that is called integral of the equation.

n ()

Let’s consider a simple ODE where is a known function, the solution is:

= (), ()

an ODE of order n has infinite solutions depending on n coefficients, the set

() = () + ,

∫

of these solutions is named general integral of the ODE.

In practical application, we’re interested in choosing one of the solutions, we can observe the

system at time t , so we know exactly y(t ) and we choose the solution that respects this condition

0 0

(it’s named particular integral).

In general, for an ODE of order n we need n initial conditions.

A Partial Differential Equation (PDE) limited to second order is:

$ $ $ $ $

•, , , ,…, , , ,…, , ,…, ,…• = 0

$

$ $

H J J H H $

H

As seen for ODE, the idea is similar, we have to include boundary conditions.

Possible approaches to solve a differential equation are analytical, qualitative and numerical ones.

Analytical solution

It’s the most powerful approach, because it allows to get the analytical expression of the solution,

but analytical solutions are available only for a limited classis of equations, probably none of the

problems in medicine and healthcare can be solved in this way.

n

If we have we have the trivial solution so y is a constant, for the other

= ()(), () = 0,

solutions we have to solve it: = ()()

‚ = ‚ () +

()

Numerical solution

It’s a flexible approach that can be always applied, but the computational time to get the solution

is high and solutions are given in terms of point without getting the analytical expression.

Numerical methods are based on the discretization of the continuous independent variable x or t in

y(x) or y(t), respectively. We consider a first order ODE:

() )

= (, ); ( =

u u

In the Euler Method, the solution at the next point is written in terms of the Taylor polynomial, cut

at the first order for approximation, so for the iterative approach, expressed in compact way we

have: )ℎ

= + ( , = 1,2, … . ℎ = ( + ℎ)

,sH , , , , u

The problem is to determine an integration step h small enough to get convergence, so we can:

1. try some values of h

2. observe that decreasing the value of h the solutions are stable

3. any value of h for which the solution does not change is suitable, the highest among them is

the most efficient from the computational viewpoint

We can also use the Runge-Kutta methods, adding additional point between x and x instead of

i i+1,

using only one slope, a set of slopes are considered to better approximate the solution; the

considered slope is then a weighted average of them.

RK numerical method are a family of iterative discrete methods utilized for the numerical

approximation of ODE’s solution, they are used to find an approximation of the function

Š that verifies the generic Cauchy problem:

(): ℝ → ℝ n () ())

= (,

‹ )

( =

u u

Generally, the problem is considered in an interval and it is considered a sample of the

[ , ]

u Œ

interval in a set of point where and then the

∆ { | = 0 … } = + ℎ ℎ = ( − )/,

, , u Πu

numerical methods gives an approximation of values.

( )

,

Given the generic Cauchi problem, with known value of and we consider a sufficiently small

u u

and positive interval and we define

ℎ ℎ

= + ( + 2 + 2 + )

,sH , H $ • •

6

= + ℎ

,sH J

In this way, is calculated as the sum of and the weighted mean of the four increments

JsH J )

= ( ,

H J J ℎ 1

= ‘ + , + ℎ’

$ J J H

2 2

ℎ 1

= ‘ + , + ℎ’

• J J $

2 2

= ( + ℎ, + ℎ)

• J J •

Then, is the increment based on the slop at the beginning of the interval using ; and are

′

H J $ •

based on the slop at the interval midpoint, using and ; is the increment based on the slop

H $ •

at the end of the interval using .

•

An ODE with n larger than one can be rewritten as an ODE with order 1, for example if we have an

ODE of order 2: n nn )

(, , , = 0

We can define two new variables, that are: ()

= ()

H ()

= ′()

$

So, the equation becomes: )

(, , , ′ = 0

H $ $

Decision, optimization approaches

Optimization means to find the maximum and minimum of a function over a given set. If the set

coincides with