Simulation of a digital blood pressure monitor and analysis of human blood pressure

Currently, there are two most common methods of measuring blood pressure - the Korotkov method (auscultation or acoustic) and oscillometric. Each of these methods provides the necessary accuracy and convenient reproducibility, however, each of them has its advantages and disadvantages.

Currently, the Korotkov method is a non–invasive blood pressure measurement method approved by the World Health Organization. The measurement is performed using a tonometer, and listening to tones from an artery clamped by a cuff is performed using a stethoscope or microphone.

The oscillometric method, like the Korotkov method, is based on recording arterial pulsations when the laminar blood flow changes to turbulent, but without listening to tones. The oscillometric method uses the recording of pressure fluctuations (oscillations) directly in the cuff. Thus, during the first Korotkov tone, the most dramatic increase in the amplitude of pulsations is observed – systolic pressure is recorded. A sharp decrease in the amplitude of pulsations indicates a change in the turbulent blood flow to laminar – diastolic pressure.

Using this example, you can study the basic principles of a digital tonometer with the oscillometric method of measuring pressure and consider the effect of human movement during the study.

Installing the necessary libraries for further modeling.

let
    installed_packages = collect(x.name for (_, x) in Pkg.dependencies() if x.is_direct_dep)
    list_packages = ["Random","Plots", "Statistics", "DSP"]
    for pack in list_packages
        pack in installed_packages || Pkg.add(pack)
    end
end

using Plots, Random, DSP, Statistics

General information about the structure of the cardiovascular system

To begin with, we will describe the physiological parameters of a person (patient) that are necessary for modeling

Pulse rate - the number of cycles of contraction (systole) and relaxation (diastole) of the heart, recorded by pulse waves in peripheral arteries per unit time
Systolic (upper) arterial pressure (SD) is the level of blood pressure at the moment of maximum contraction of the heart, characterizes the state of the myocardium of the left ventricle.
Diastolic (lower) blood pressure (DB) is the blood pressure level at the moment of maximum relaxation of the heart, characterizes the degree of tone of the arterial walls.

# Human physiological parameters (paceinta)
age = 20            # Age, year
weight = 80         # Weight, kg
height = 180        # Height, cm

pulse_rate = 75;             # Pulse rate, beats/min.
systolic_bp = 120.0;         # Systolic blood pressure, mmHg
diastolic_bp = 80.0;         # Diastolic pressure, mmHg

Blood pressure is measured in millimeters of mercury, abbreviated as mmHg. A blood pressure value of 120/80 means that systolic pressure is 120 mmHg, and diastolic blood pressure is 80 mmHg.

The difference between systolic and diastolic pressure values is called pulse pressure (PD). It shows how much systolic pressure exceeds diastolic pressure, which is necessary to open the semilunar aortic valve during systole. Normally, the pulse pressure is 35-55 mmHg.

Blood pressure (BP) is measured in millimeters of mercury (mmHg).
Writing a blood pressure value of 120/80 mmHg means:

Systolic blood pressure (DM) = 120 mmHg.
Diastolic blood pressure (DD) = 80 mmHg.

Pulse pressure (PD) is the difference between systolic and diastolic pressure values.:

pulse_bp = systolic_bp - diastolic_bp;      # Pulse pressure, mmHg

To determine the average arterial pressure (MPP), which expresses the energy of continuous blood flow, the following Hickam formula can be used:

mean_bp = (systolic_bp + 2*diastolic_bp)/3;   # Average blood pressure, mmHg

To assess the functional state of the cardiovascular system, the minute volume of the heart (MO) is calculated and compared with the proper value (DMO).

where 2.2 is the cardiac index, measured in liters, and PT is the body surface, calculated using a nomogram or various formulas. In this example, the Dubois formula will be used.:

And the stroke volume (UO) of the heart is calculated by the expression:

where B is the age, g.

Let's introduce the notation in English

PT (body surface) → BSA (Body Surface Area)
UO (impact volume) → SV (Stroke Volume)
MO (minute volume) → CO (Cardiac Output)
DMO (due minute volume) → RCO (Required Cardiac Output)

BSA =  0.007184 * (weight^0.425) * (height^0.725);                  # Body surface, m^2
SV = 101 + 0.5 * systolic_bp - 1.08 * diastolic_bp - 0.6 * age;     # Impact volume, ml
CO = SV * pulse_rate;                                               # Minute volume, ml/min
RCO = 2.2 * BSA * 1000;                                             # Required minute volume, ml/min

We will display the obtained results of calculating the parameters of the cardiovascular system (CVS).

println("The results of calculating the parameters of the cardiovascular system (CVS):")
println("----------------------------------------------------------------")
println("1. Pulse pressure (PP): ", round(pulse_bp; digits=1), " mmHg")
println("2. Average arterial pressure (MPP): ", round(mean_bp; digits=1), " mmHg")
println("3. Body Surface Area (PT): ", round(BSA; digits=4), " m2")
println("4. Impact volume (UO): ", round(SV; digits=1), " ml")
println("5. Minute Heart Volume (MO): ", round(CO; digits=1), " ml/min")
println("6. Proper Minute Volume (DMO): ", round(RCO; digits=1), " ml/min")

Результаты расчета параметров сердечно сосудистой системы (ССС):
----------------------------------------------------------------
1. Пульсовое давление (ПП): 40.0 мм.рт.ст
2. Среднее артериально давление (СрАД): 93.3 мм.рт.ст
3. Площадь поверхности тела (ПТ): 1.9964 м²
4. Ударный объем (УО): 62.6 мл
5. Минутный объем сердца (МО): 4695.0 мл/мин
6. Должный минутный объем (ДМО): 4392.1 мл/мин

Simulation of the operation of a digital tonometer with an oscillometric measurement method

A digital tonometer with oscillometric measurement method determines blood pressure by analyzing fluctuations (oscillations) of air pressure that occur in the cuff when a pulse wave is transmitted from the artery. Let's analyze the stages of his work sequentially.:

Increasing the pressure in the cuff
Slow release of pressure in the cuff
Registration of fluctuations
Processing of measurement results
Displaying the results

The preparatory stage of modeling consists in setting the signal parameters necessary for the subsequent generation of physiological data. This stage determines the key characteristics of the time axis on which all signals will be built.

# Signal Parameters
fs = 200            # Sampling rate (200 Hz)
total = 60.0        # Total measurement duration (60 seconds)

# Timeline
t = 0:1/fs:total - 1/fs;

The first stage is to simulate the pressure in the cuff of the tonometer. The model uses three key parameters to generate a realistic pressure curve.:

Duration of pumping: 10.0 seconds
(optimal time for comfortable cuff filling without excessive load on the arteries)
Plateau duration: 1.0 sec
(the period of pressure stabilization for complete compression of the artery)
Peak pressure: 160.0 mmHg
(standard value exceeding typical systolic pressure by 30-40 mmHg)

# Phases of pressure measurement in the tonometer cuff
time_pump = 10.0                # Duration of cuff inflation, seconds
time_plato = 1.0                # Plateau duration after pumping, seconds
peak_pressure = 160.0           # Maximum pressure in the cuff, mmHg


# Simulation of pressure in the tonometer cuff
function generate_pressure(time)
    pressure = zeros(length(time))       # Creating an array of zeros of the same length as the time vector
    
    # Pumping pressure into the cuff
    pump_phase = time .<= time_pump
    pressure[pump_phase] = peak_pressure * (time[pump_phase] / time_pump).^0.85

    # The plateau
    plato_start = time_pump
    plato_end = time_pump + time_plato
    plato_phase = (time .> plato_start) .& (time .<= plato_end)
    pressure[plato_phase] .= peak_pressure

    # Cuff pressure relief
    deflate_phase = (time.> plato_end)
    deflate_time = time[deflate_phase] .- plato_end
    tau = 7.0 
    pressure[deflate_phase] = peak_pressure .* exp.(-deflate_time ./ tau)
    
    return pressure, deflate_phase, plato_end
end;

The second stage consists in modeling the oscillations in the cuff of the tonometer. The generart_oscillation function implements fluctuations taking into account the following aspects:

The systolic phase corresponds to a direct pulse beat.
The decrotic recession simulates a reflected wave from a closed aortic valve.
Diastolic fluctuations reflect residual pulsation in blood vessels.
Gaussian modulation simulates a real phenomenon: the maximum amplitude of oscillations is observed when the pressure in the cuff is equal to the average arterial pressure.

function generate_oscillations(pressure, time, pulse_rate, systolic_bp, diastolic_bp, mean_bp, inflation_duration; movement=false)
    pulse_freq = pulse_rate / 60.0
    T = 1 / pulse_freq
    oscillations = zeros(length(time))
    movement_detected = false
    movement_start_time = 0
    movement_end_time = 0
    movement_start = 140.0          # Initial pressure of motion artifacts, mmHg
    movement_end = 110.0            # The final pressure of motion artifacts, mmHg
    
    for i in eachindex(time)
        t_val = time[i]
        p = pressure[i]
        t_mod = mod(t_val, T)
        
        # The phase of the cardiac cycle
        systolic_phase = 0.15 * T       # Systole, sec
        dicrotic_phase = 0.18 * T       # Decrotic notch, sec
        diastolic_phase = 0.45 * T      # Diastole, sec
        
        wave = 0.0
        if t_mod < systolic_phase
            phase_ratio = t_mod / systolic_phase
            wave = 1.6 * phase_ratio * exp(1.8 * (1 - phase_ratio))
        elseif t_mod < dicrotic_phase
            phase_ratio = (t_mod - systolic_phase) / (dicrotic_phase - systolic_phase)
            wave = 1.0 - 2.2 * phase_ratio
        elseif t_mod < systolic_phase + diastolic_phase
            phase_ratio = (t_mod - dicrotic_phase) / (diastolic_phase - (dicrotic_phase - systolic_phase))
            wave = 0.3 * exp(-3.5 * phase_ratio) * sin(2*π * 1.2 * phase_ratio)
        end
        
        # Amplitude modulation
        spread = (systolic_bp - diastolic_bp) / 2.5
        amp_mod = exp(-0.5 * ((p - mean_bp) / spread)^2)
        
        if t_val < inflation_duration
            amp_mod = 0.0
        elseif p < 60
            amp_mod = 0.0
        elseif p > systolic_bp + 20 || p < diastolic_bp - 15
            amp_mod *= 0.02
        elseif p < 80
            amp_mod *= (p - 60) / 20
        end
        
        # Adding basic noise
        noise = 0.03 * randn()
        oscillations[i] = 5.2 * amp_mod * (wave + noise)
        
        # Generating motion artifacts
        if movement
            if t_val > inflation_duration && p <= movement_start && p >= movement_end
                if !movement_detected
                    movement_detected = true
                    movement_start_time = t_val
                end
                oscillations[i] += 3.5 * sin(2π * 0.8 * t_val)
                if rand() < 0.15
                    oscillations[i] += 8.0 * randn()
                end
            else
                if movement_detected
                    movement_detected = false
                    movement_end_time = t_val
                end
            end
        end
    end
    
    # If the movement has not ended
    if movement && movement_detected
        movement_end_time = time[end]
    end
    
    return oscillations, movement_start_time, movement_end_time
end;

The third stage consists in analyzing the results obtained and detecting DM, DD and SrAD using generally accepted methods.

function analysis(raw_osc, deflate_phase, pressure, time, fs)
    # Signal filtering
    ttime = 0.15                        # Filter reaction time
    n = Int(round(ttime * fs))          # Filter window size
    b = ones(n)/n                       # FIR filter coefficients
    filt_osc = filt(b, [1], raw_osc)    # Applying the filter
    
    # Detection of the pressure deflation phase in the cuff
    deflate_osc = filt_osc[deflate_phase]
    deflate_pres = pressure[deflate_phase]
    deflate_times = time[deflate_phase]
    
    # Building an envelope
    abs_osc = abs.(deflate_osc)
    window_size = Int(round(0.45 * fs))
    window = ones(window_size)/window_size
    envelope = filt(window, [1], abs_osc)
    
    # We find the maximum amplitude
    max_amp, max_idx = findmax(envelope)
    pres_max = deflate_pres[max_idx]
    
    # Systolic blood pressure
    systolic_lim = 0.15 * max_amp
    systolic_region = findfirst(x -> x >= systolic_lim, envelope[1:max_idx])
    systolic_idx = something(systolic_region, 1)
    systolic = deflate_pres[systolic_idx]
    
    # Diastolic pressure
    diastolic_lim = 0.55 * max_amp
    diastolic_region = findlast(x -> x >= diastolic_lim, envelope[max_idx:end])
    diastolic_idx = diastolic_region === nothing ? length(envelope) : max_idx - 1 + diastolic_region
    diastolic = deflate_pres[diastolic_idx]
    
    return deflate_osc, deflate_times, envelope, systolic, diastolic, pres_max
end;

The first scenario.

In this scenario, we will consider the case when the patient adheres to the rules for measuring blood pressure.

# First scenario: no motion artifacts
pressure, deflate_phase, plato_end = generate_pressure(t)
time_period = time_pump + time_plato
osc, _, _ = generate_oscillations(pressure, t, pulse_rate, systolic_bp, diastolic_bp, mean_bp, time_period, movement=false)
deflate_osc, deflate_times, envelope, systolic_est, diastolic_est, pres_max = analysis(osc, deflate_phase, pressure, t, fs)
pres_osc = pressure .+ osc

# Creating a graph with two subgraphs
plt_1 = plot(layout=(2,1), size=(900, 700), legend=:topright, margin=40*Plots.px)

# Upper graph: Pressure in the tonometer cuff
plot!(plt_1[1], t, pres_osc, 
    title = "Pressure in the tonometer cuff (without artifacts)",
    xlabel = "Time, from",
    ylabel = "Pressure, mmHg",
    label = "Cuff pressure",
    linewidth = 2,
    color = :blue,
    ylims = (0, 200),
    xlims = (0, total),
    grid = true)

# Key pressure points
hline!(plt_1[1], [systolic_est], linestyle=:dash, color=:red, 
       label="Systolic (ism: $(round(systolic_est, digits=1))")
hline!(plt_1[1], [diastolic_est], linestyle=:dash, color=:purple, 
       label="Diastolic (ism: $(round(diastolic_est, digits=1))")
hline!(plt_1[1], [pres_max], linestyle=:dash, color=:brown, 
       label="Average (ism: $(round(pres_max, digits=1))")

# Real values
hline!(plt_1[1], [systolic_bp], linestyle=:dot, color=:red, label="Systolic (real: $systolic_bp)")
hline!(plt_1[1], [diastolic_bp], linestyle=:dot, color=:purple, label="Diastolic (real: $diastolic_bp)")

# Lower graph: Pressure fluctuations
plot!(plt_1[2], t, osc, 
    title = "Pressure fluctuations",
    xlabel = "Time, from",
    ylabel = "Amplitude, mmHg",
    label = "The oscillations",
    linewidth = 1.5,
    color = :blue,
    xlims = (0, total),
    grid = true)

# The envelope is only in the deflation phase
plot!(plt_1[2], deflate_times, envelope, 
    label = "The oscillation envelope",
    linewidth = 2.0,
    color = :red)

# Definition thresholds
systolic_lim = 0.15 * maximum(envelope)
diastolic_lim = 0.55 * maximum(envelope)
hline!(plt_1[2], [systolic_lim], linestyle=:dash, color=:green, label="Systolic threshold (15%)")
hline!(plt_1[2], [diastolic_lim], linestyle=:dash, color=:orange, label="Diastolic threshold (55%)")
hline!(plt_1[2], [maximum(envelope)], linestyle=:dash, color=:purple, label="Maximum amplitude")

# Output of the results of the first scenario
println("Results of the analysis of the first scenario (lack of movement):")
println("Systolic blood pressure: $(round(systolic_est, digits=1)) mmHg (real: $systolic_bp)")
println("Diastolic pressure: $(round(diastolic_est, digits=1)) mmHg (real: $diastolic_bp)")
println("Average blood pressure: $(round(pres_max, digits=1)) mmHg")
println("Pulse rate: $pulse_rate beats/min")

# Graph displays
display(plt_1)
savefig(plt_1, "tonometer_no_movement.png")

Результаты анализа первого сценария (отсутствие движения):
Систолическое давление: 120.4 мм рт.ст. (реальное: 120.0)
Диастолическое давление: 73.2 мм рт.ст. (реальное: 80.0)
Среднее артериальное давление: 92.2 мм рт.ст.
Пульс: 75 уд/мин

"/user/Demo_public/biomedical/modeling_blood_pressure_monitor/tonometer_no_movement.png"

From the simulation results of the first scenario, it can be seen that the resulting model detects fluctuations in human blood pressure with an accuracy of at least 5%.

The second scenario.

In this scenario, we will consider the case when the patient ignores the rules for measuring blood pressure.

# The second scenario: with motion artifacts
osc_movement, movement_start_time, movement_end_time = generate_oscillations(pressure, t, pulse_rate, systolic_bp, diastolic_bp, mean_bp, time_period, movement=true)
deflate_osc_movement, deflate_times_movement, envelope_movement, systolic_est_movement, diastolic_est_movement, max_pres_movement = analysis(osc_movement, deflate_phase, pressure, t, fs)
pressure_movement = pressure .+ osc_movement

# Creating a graph (two subgraphs)
plt_2 = plot(layout=(2,1), size=(900, 700), legend=:topright, margin=40*Plots.px)

# Upper graph: Tonometer cuff pressure with artifacts
plot!(plt_2[1], t, pressure_movement, 
    title = "Pressure in the tonometer cuff (with motion artifacts)",
    xlabel = "Time, from",
    ylabel = "Pressure, mmHg",
    label = "Cuff pressure",
    linewidth = 2,
    color = :blue,
    ylims = (0, 200),
    xlims = (0, total),
    grid = true)

# Highlighting the area of motion artifacts
if !isnan(movement_start_time) && !isnan(movement_end_time)
    vspan!(plt_2[1], [movement_start_time, movement_end_time], alpha=0.25, color=:magenta, label="Movement artifacts")
end

# Key pressure points
hline!(plt_2[1], [systolic_est_movement], linestyle=:dash, color=:red, label="Systolic (ism: $(round(systolic_est_movement, digits=1))")
hline!(plt_2[1], [diastolic_est_movement], linestyle=:dash, color=:purple,label="Diastolic (ism: $(round(diastolic_est_movement, digits=1))")
hline!(plt_2[1], [max_pres_movement], linestyle=:dash, color=:brown, label="Average (ism: $(round(max_pres_movement, digits=1))")

# Real pressure values
hline!(plt_2[1], [systolic_bp], linestyle=:dot, color=:red, label="Systolic (real: $systolic_bp)")
hline!(plt_2[1], [diastolic_bp], linestyle=:dot, color=:purple, label="Diastolic (real: $diastolic_bp)")

# Bottom chart: Oscillations with motion artifacts
plot!(plt_2[2], t, osc_movement, 
    title = "Pressure fluctuations with artifacts",
    xlabel = "Time, from",
    ylabel = "Amplitude, mmHg",
    label = "The oscillations",
    linewidth = 1.5,
    color = :blue,
    xlims = (0, total),
    grid = true)

# The envelope
plot!(plt_2[2], deflate_times_movement, envelope_movement, 
    label = "The envelope",
    linewidth = 2.0,
    color = :red)

# Definition thresholds
systolic_lim_mov = 0.15 * maximum(envelope_movement)
diastolic_lim_mov = 0.55 * maximum(envelope_movement)
hline!(plt_2[2], [systolic_lim_mov], linestyle=:dash, color=:green, label="Systolic threshold (15%)")
hline!(plt_2[2], [diastolic_lim_mov], linestyle=:dash, color=:orange, label="Diastolic threshold (55%)")
hline!(plt_2[2], [maximum(envelope_movement)], linestyle=:dash, color=:purple, label="Maximum amplitude")

# The area of motion artifacts
if !isnan(movement_start_time) && !isnan(movement_end_time)
    vspan!(plt_2[2], [movement_start_time, movement_end_time], 
           alpha=0.15, color=:magenta, label="Movement artifacts")
end

# Output of results for the second scenario
println("\Results of the analysis of the second scenario (presence of movement):")
println("Systolic pressure: $(round(systolic_est_movement, digits=1)) mmHg (real: $systolic_bp)")
println("Diastolic pressure: $(round(diastolic_est_movement, digits=1)) mmHg (real: $diastolic_bp)")
println("Average blood pressure: $(round(max_pres_movement, digits=1)) mmHg")
println("Pulse rate: $pulse_rate beats/min")

# Graph display
display(plt_2)
savefig(plt_2, "tonometer_with_movement.png")

Результаты анализа второго сценария (наличие движения):
Систолическое давление: 137.0 мм рт.ст. (реальное: 120.0)
Диастолическое давление: 81.3 мм рт.ст. (реальное: 80.0)
Среднее артериальное давление: 118.4 мм рт.ст.
Пульс: 75 уд/мин

"/user/Demo_public/biomedical/modeling_blood_pressure_monitor/tonometer_with_movement.png"

According to the results of the second scenario, it is clear that the determination of systolic pressure is not correct. Erroneous detection occurs due to artifacts caused by human movement during measurement.

Visualization of scenarios for measuring human pressure in the absence and presence of artifacts caused by human movement.

# Visualization of the results of modeling the first and second scenarios
# Output of the results of the first scenario
println("Results of the analysis of the first scenario (lack of movement):")
println("Systolic blood pressure: $(round(systolic_est, digits=1)) mmHg (real: $systolic_bp)")
println("Diastolic pressure: $(round(diastolic_est, digits=1)) mmHg (real: $diastolic_bp)")
println("Average blood pressure: $(round(pres_max, digits=1)) mmHg")
println("Pulse rate: $pulse_rate beats/min")

# Graph displays
display(plt_1)

# Output of results for the second scenario
println("\Results of the analysis of the second scenario (presence of movement):")
println("Systolic pressure: $(round(systolic_est_movement, digits=1)) mmHg (real: $systolic_bp)")
println("Diastolic pressure: $(round(diastolic_est_movement, digits=1)) mmHg (real: $diastolic_bp)")
println("Average blood pressure: $(round(max_pres_movement, digits=1)) mmHg")
println("Pulse rate: $pulse_rate beats/min")

# Graph display
display(plt_2)

Результаты анализа первого сценария (отсутствие движения):
Систолическое давление: 120.4 мм рт.ст. (реальное: 120.0)
Диастолическое давление: 73.2 мм рт.ст. (реальное: 80.0)
Среднее артериальное давление: 92.2 мм рт.ст.
Пульс: 75 уд/мин

Результаты анализа второго сценария (наличие движения):
Систолическое давление: 136.3 мм рт.ст. (реальное: 120.0)
Диастолическое давление: 81.1 мм рт.ст. (реальное: 80.0)
Среднее артериальное давление: 129.0 мм рт.ст.
Пульс: 75 уд/мин

Conclusion

In this example, a digital tonometer model with an oscillometric pressure measurement method was considered. Two scenarios were also considered:

cuff pressure without artefacts of human (patient) movement;
pressure in the cuff with artifacts of human (patient) movement.

Analyzing the results obtained, it can be said that in order to correctly determine a person's systolic and diastolic blood pressure, a number of rules must be followed.:

the measurement should be carried out in a comfortable, calm environment, the room should be at room temperature;
Blood pressure can only be measured after at least a five-minute rest period.;
it should be remembered that the shoulder should not be squeezed by clothes, especially since it is incorrect to measure blood pressure through clothes;
Do not move or talk during the measurement.

During the initial measurement, blood pressure should be determined on both hands and then blood pressure should be measured on the hand where the pressure was higher. (The difference in blood pressure on the hands up to 10-15 mmHg is normal.)