A sensitivity-optimized population and individualized respiratory process modeling method

Abstract

The application discloses a kind of group and individualized respiratory process modeling method based on sensitivity optimization, belong to biomedical engineering and physiology technical field, including the following steps: S1, data preparation;S2, model construction and parameter optimization;S3, model verification.The application adopts the above-mentioned group and individualized respiratory process modeling method based on sensitivity optimization, which realizes the accurate simulation of blood oxygen saturation, respiratory rate and tidal volume, captures the hypoxic response characteristics, avoids the influence of low-quality signals, and makes the daily difference simulation have higher fitting accuracy.

Description

TECHNICAL FIELD

[0001] The present application belongs to the technical field of biomedical engineering and physiology, and particularly relates to a population and individualized respiratory process modeling method based on sensitivity optimization. BACKGROUND

[0002] Respiration is one of the important physiological functions of the human body, and ensuring normal respiration of the human body is the most basic requirement for maintaining life. When a person on the plain is in a highland environment, respiration will be limited by the decrease of oxygen concentration, which will cause acute high altitude reaction and even induce high altitude mountain sickness. Therefore, it is very valuable to construct a physiological model that can predict indicators such as blood oxygen saturation (SpO2), respiratory rate, and tidal volume for high altitude disease prevention, low oxygen preparation, etc.

[0003] At present, many respiratory models based on physiological principles have been established. The main dynamic characteristics of respiration can be described in the following three aspects. First, lung gas exchange and blood perfusion. Second, the central nervous system (CNS) controls the respiratory process. The CNS model of the respiratory process is constructed based on neuron electrical signal conduction, which accurately describes the working mechanism of the respiratory mode and the brainstem. Third, the cardiovascular system regulates respiration. The oxygen dissociation curve fitting equation is used to describe the change of blood hemoglobin (Hb) oxygen affinity with oxygen partial pressure.

[0004] Population and individualized model construction requires an effective parameter estimation scheme. Variational inference (VI) is widely used to evaluate complex injuries, and the key is to select a distribution close to the true posterior distribution. The prior art optimizes the objective through collective variational inference (C-VI) to coordinate low-quality and heterogeneous data, but the signal quality measurement standard is not established.

[0005] Because the existing model faces the challenges of limited sample size and uneven signal quality. Therefore, there is an urgent need for a method to construct an individualized physiological model. SUMMARY

[0006] The purpose of the present application is to provide a population and individualized respiratory process modeling method based on sensitivity optimization, which realizes accurate simulation of blood oxygen saturation, respiratory rate and tidal volume, captures hypoxia reaction characteristics, avoids the influence of low-quality signals, and makes the daily difference simulation have higher fitting accuracy.

[0007] To achieve the above purpose, the present application provides a population and individualized respiratory process modeling method based on sensitivity optimization, comprising the following steps:

[0008] S1, collecting data and preprocessing the data;

[0009] S2, initialize the population and individualized human respiratory process physiological model in the first module, bring the data in S1 into the second module as real data, optimize the parameters based on the identifiable integrated variational inference and hypoxia perception maximum likelihood estimation, use the sensitivity-based swarm intelligence genetic algorithm in the third module to optimize the second module, and obtain the optimized population and individualized parameters;

[0010] S3, use the test data set to verify the prediction accuracy of the population and individualized human respiratory process physiological model constructed in S2.

[0011] A population and individualized respiratory process modeling system based on sensitivity optimization, comprising: a preprocessing module for data acquisition and preprocessing; a parameter processing module for initializing the model and optimizing the data to obtain optimized population and individualized parameters.

[0012] Preferably, the parameter processing module comprises: a first module for constructing a population and individualized human respiratory process physiological model; a second module using identifiable integrated variational inference and hypoxia perception maximum likelihood estimation for parameter estimation; and a third module using a sensitivity-based swarm intelligence genetic algorithm to optimize the objective function in the second module.

[0013] Preferably, the first module comprises a respiratory system model, a central nervous system model, and a cardiovascular system model.

[0014] Preferably, the central nervous system model simplifies the transmission of neural signals in the brainstem and designs a segmented proportional-integral controller to describe the control of respiratory frequency and tidal volume changes through the following equations:

[0015] R r = R r0 +k pr (S r -S e )+k ir ∫(S r -S e )dt (5)

[0016] ΔV t = ΔV t0 +k pv (S r -S e )+k iv ∫(S r -S e )dt (6)

[0017]

[0018] wherein, R r represents the respiratory frequency, and ΔVt represents tidal volume variation, S e represents Sp02monitored by chemical sensor, S r represents blood oxygen saturation reference value, (S r -S e ) represents deviation value of feedback signal, k pr represents proportional controller parameter of respiratory frequency, k pv represents integral control parameter of respiratory frequency, k ir represents proportional control parameter of respiratory amplitude, k iv represents integral controller parameter of respiratory amplitude, R r0 represents initial value of respiratory frequency, ΔV t0 represents initial value of tidal volume variation.

[0019] Preferably, the human respiratory process physiological model is described by the following equations:

[0020] x = [P A , S a , S v , R r , ΔV t ] (9)

[0021] u = [P m , H r ] (10)

[0022] y = S a (11)

[0023] θ = [k1, k2, k3, k4, w, V r , V t , S r , k pr1 , k pr2 ,

[0024] k ir1 , k ir2 , k pv1 , k pv2 , k iv1 , k iv2 , P 50 , r] (12)

[0025] wherein x represents model state, u represents input data, y represents output data, θ represents key parameter vector, P A represents alveolar oxygen partial pressure, S a represents blood oxygen saturation, k1 represents diffusion ability of oxygen from alveoli to arterial blood, k2 represents diffusion ability of oxygen from arterial blood to tissue, k3 represents diffusion ability of oxygen from tissue to venous blood, k4 represents diffusion ability of oxygen from venous blood to alveoli, S vdenotes venous oxygen saturation, R r denotes respiratory rate, P m denotes inhaled gas oxygen partial pressure, H r denotes heart rate, w denotes the effect of metabolism on venous oxygen saturation, V r denotes functional residual capacity, V t denotes baseline tidal volume, P 50 denotes oxygen partial pressure at which blood oxygen saturation is 50%, r denotes the steepness of the oxygen dissociation curve.

[0026] Preferably, the diurnal variation parameter θ v = [k1, k2, P 50 , r] varies over time.

[0027] Preferably, in the second module, the optimization objective based on identifiable integrated variational inference is represented as:

[0028] L(v) = log p(y) - D KL (v)

[0029] = E q [log p(y|θ, u) + log p(θ|φ) + log p(φ) - log q(φ, θ|v)] (16)

[0030] wherein log p(y|θ, u) is represented as:

[0031]

[0032] l ij denotes signal length, denotes simulated output data, denotes real output data, u denotes input data;

[0033] log q(φ, θ|v) is represented as:

[0034]

[0035] ∈ k denotes the error of the kth group, v σ denotes the precision parameter of the variational distribution of the kth group;

[0036] log p(θ|φ) is represented as:

[0037]

[0038] ||.|| * denotes the kernel norm of the matrix, λ denotes a hyperparameter;

[0039] log p(φ) is represented as:

[0040] log p(φ) = -λ||φ ∑ || * (23)

[0041] φ μ denotes the mean of the individualized parameters, φ Σ denotes the correlation between the individualized parameters, φ denotes the variational distribution.

[0042] Preferably, in the second module, an anoxic-aware maximum likelihood estimation is added in the parameter estimation, which is calculated as follows:

[0043]

[0044] where ∑ w denotes the diagonal weight matrix, θ i denotes the individualized parameters of the i-th experimental subject, u i denotes the input data of the experimental subject i.

[0045] Preferably, in the third module, the sensitivity-based swarm intelligence genetic algorithm includes the following steps:

[0046] Step 1, initialize the population, randomly generate multiple individuals containing complete parameter vectors;

[0047] Step 2, according to the order of magnitude of the parameter sensitivity, divide the parameters into a high-sensitivity group g high and a low-sensitivity group g low ;

[0048] Step 3, calculate the fitness of each individual using the identifiable-based integrated variational inference objective function to evaluate its pros and cons;

[0049] Step 4, apply the genetic algorithm to the high-sensitivity group g high , generate new individuals by selecting individuals with high fitness for crossover and mutation operations;

[0050] Step 5, apply the genetic algorithm to the low-sensitivity group g low , generate new individuals by selecting individuals with high fitness for crossover and mutation operations;

[0051] Step 6, use the set Ω + to select high-sensitivity parameters, and extract the gene group of high-sensitivity parameters from the chromosome for optimization, and then restore the gene group of high-sensitivity parameters to the entire chromosome, use the set Ω - to select low-sensitivity parameters, and extract the gene group of low-sensitivity parameters from the chromosome for optimization, and then restore the gene group of low-sensitivity parameters to the entire chromosome;

[0052] Step 7, repeat the above steps until a preset number of iterations is reached or a termination condition is met.

[0053] Therefore, the present application adopts the above-mentioned group and individualized respiratory process modeling method based on sensitivity optimization, compared with the prior art, the present application has the following remarkable beneficial effects:

[0054] (1) The present application realizes the simulation of blood oxygen saturation, respiratory rate and tidal volume and other indicators by including the control of the nervous system on the respiratory system in modeling, and captures the characteristics of acute and chronic hypoxia reactions;

[0055] (2) The present application designs an IBC-VI optimization target based on the recognition, innovatively measures the signal quality, avoids the negative impact of low-quality signals on distribution estimation, and improves the adaptability of the group and individual multi-level evaluation model;

[0056] (3) The present application proposes an optimization scheme based on sensitivity analysis (SA) and a group intelligence genetic algorithm based on sensitivity (SGBW-GA), which calculates the sensitivity of the model parameters and optimizes the parameters with large differences respectively, effectively improving the optimization precision;

[0057] (4) The present application proposes a multi-parameter estimation optimization scheme to improve data utilization and solve the limitations of small sample learning;

[0058] (5) The diurnal difference simulation of the present application has higher fitting precision, and can more clearly understand the complex physiological process of the human body and analyze the similarities and differences between people.

[0059] The technical solutions of the present application will be further described in detail below with the help of the drawings and examples. BRIEF DESCRIPTION OF DRAWINGS

[0060] Figure 1 The SpO2 data graph of the present application is based on a sensitivity optimization group and individualized respiratory process modeling method, each subgraph includes three groups, which respectively represent the measured value, the result of the SBGW-GA method, the result of the basic GA method and the result of the diurnal difference parameter, a represents the data of the first day, the MSE of the SBGW-GA method is 33.7, the MSE of the basic GA method is 21.82, and the MSE of the diurnal difference parameter is 20.11; b represents the data of the second day; c represents the data of the fifth day; d represents the data of the sixth day; e represents the data of the ninth day; f represents the data of the tenth day;

[0061] Figure 2The respiratory frequency curve chart of the population and individualized respiratory process modeling method based on sensitivity optimization, each subgraph includes three groups, the three groups from top to bottom represent the real curve, the SBGW-GA simulation curve and the basic GA simulation curve, a represents the change of the respiratory frequency of the three different groups with time at the second day; b represents the change of the respiratory frequency of the three different groups with time at the sixth day; c represents the change of the respiratory frequency of the three different groups with time at the tenth day;

[0062] Figure 3 The tidal volume curve chart of the population and individualized respiratory process modeling method based on sensitivity optimization, each subgraph includes three groups, the three groups from top to bottom represent the basic GA simulation curve, the diurnal difference parameter result simulation curve and the SBGW-GA simulation curve, a represents the change of the tidal volume of the three different groups with time at the second day; b represents the change of the tidal volume of the three different groups with time at the sixth day; c represents the change of the tidal volume of the three different groups with time at the tenth day;

[0063] Figure 4 The box chart of the individualized parameter and distribution sampling of the population and individualized respiratory process modeling method based on sensitivity optimization, each subgraph includes three different groups, the three groups from left to right are the individualized parameter θ i , the sampling data of φ obtained by using the IBC-VI optimization target and the sampling data of φ obtained by using the C-VI optimization target, a represents the difference of the parameter k1 between different groups; b represents the difference of the parameter k2 between different groups; c represents the difference of the parameter P 50 between different groups; d represents the difference of the parameter r between different groups, "+" represents an abnormal value;

[0064] Figure 5 The physiological parameter simulation curve chart of the population and individualized respiratory process modeling method based on sensitivity optimization, wherein a shows the change of the blood oxygen saturation at different time points, the mean value of IBC-VI gradually rises; b shows the change of the tidal volume with time, the mean value of the tidal volume gradually rises with the lapse of time, and the mean value of IBC-VI is always lower than that of C-VI; c shows the change of the respiratory frequency with time, the mean value of the respiratory frequency also gradually rises with the lapse of time, and the mean value of IBC-VI is always higher than that of C-VI;

[0065] Figure 6 The diurnal difference parameter θ vThe change curve diagram shows that each dotted line represents the data curve of an experimental object, and the thick solid lines represent the data curves of experimental objects 4 and 6 respectively. a shows the change of parameter k1 for each experimental object over a period of time; b shows the change of parameter k2 for each experimental object over a period of time; c shows the change of parameter P for each experimental object over a period of time. 50 ; d shows the change of parameter r for each experimental object over a period of time. DETAILED DESCRIPTION

[0066] In order to make the purpose, technical solutions and advantages of the embodiments of the present invention clearer, the technical solutions in the embodiments of the present invention will be clearly and completely described below in conjunction with the drawings in the embodiments of the present invention. Obviously, the described embodiments are part of the embodiments of the present invention, not all of the embodiments. All other embodiments obtained by ordinary technicians in this field without making creative work are within the scope of protection of the present invention. Unless otherwise defined, the technical terms or scientific terms used in the present invention should be the common meanings understood by people with ordinary skills in the field to which the present invention belongs.

[0067] Example 1

[0068] like Figure 1 As shown, a group and individual respiratory process modeling method based on sensitivity optimization of the present invention includes the following steps:

[0069] S1. Use the preprocessing module to collect data and preprocess the data;

[0070] In this embodiment, intermittent hypoxia training (IHT) data is used. This data set contains P m 、S a and H r Data, where P m Indicates the partial pressure of oxygen in the inspired gas, S a Indicates blood oxygen saturation, H r Indicates heart rate. The experimental subjects were 19 healthy subjects aged between 18 and 30 years old, and each subject underwent IHT for 9-14 days. Daily IHT adjusted the concentration of inhaled oxygen so that the subjects were in a cyclical state of 5 minutes of hypoxia (11%-14% O2) and 3 minutes of hyperoxia (35%-38% O2), with a total duration of 1 hour. In addition, a hypoxia tolerance test was conducted on the first and last day of IHT to reflect the changes in the subjects' hypoxia adaptability. The change in hypoxia tolerance index is expressed as ΔHTI.

[0071] The preprocessing module discretizes the differential equation using the Euler method with a sampling time of T = 0.01s and preprocesses the data set using piecewise cubic Hermite interpolation.

[0072] S2. The parameter processing module includes: the first module is used to construct group and individualized physiological models of the human respiratory process; the second module uses identifiability-based integrated variational inference and hypoxia perception maximum likelihood estimation to optimize parameters; the third module uses a sensitivity-based swarm intelligence genetic algorithm to optimize the second module.

[0073] The first module is used to construct group and individualized physiological models of the human respiratory process. The model consists of three parts: the respiratory system model, the central nervous system model, and the cardiovascular system model.

[0074] The respiratory system mainly realizes the pulmonary gas exchange of the human body. Pulmonary gas exchange refers to the process of renewing oxygen in the human body and expelling carbon dioxide. During the inhalation process, oxygen diffuses from the outside into the alveoli, and then diffuses from the alveoli into the arteries. However, in a hypoxic environment, the diffusion of gas is limited by air pressure, which leads to a decrease in the oxygen content in human tissues. In order to reduce this negative impact, the body will spontaneously produce an acute hypoxic ventilation response. This response refers to the increase in respiratory ventilation in the first few minutes when a person is exposed to a hypoxic environment. After a few tens of minutes after exposure, this response will also weaken. In order to describe the oxygen diffusion process and consider that the lung volume is controlled by the nervous system, a respiratory system model is constructed, as follows:

[0075]

[0076]

[0077] Among them, P A Indicates alveolar oxygen partial pressure in mmHg. P A Derivative, R r (V t +ΔV t ) represents respiratory ventilation, the unit is L / min, V L Indicates lung volume in L, P m Indicates the oxygen partial pressure of the inhaled gas in mmHg, S a Indicates blood oxygen saturation (SpO2), S a Derivative, S v Indicates venous oxygen saturation (SvO2), S vderivative, w represents the effect of metabolism on SvO2, with units of 1 / s, k1 represents the diffusing capacity of oxygen from alveoli to arterial blood, k2 represents the diffusing capacity of oxygen from arterial blood to tissue, k3 represents the diffusing capacity of oxygen from tissue to venous blood, k4 represents the diffusing capacity of oxygen from venous blood to alveoli, k1 has units of mmHg / s, k2, k3, and k4 have units of 1 / s, f d represents the oxygen dissociation curve function.

[0078] The central nervous system's role in respiration is mainly to control pulmonary gas exchange. Central pattern generators in the brainstem control the respiratory rhythm by receiving deviation signals in mechanical and chemical feedback pathways. In mechanical feedback, the central nervous system controls normal expiration and inspiration based on the changes in lung volume detected by pulmonary mechanoreceptors. In the respiratory system model, the asymmetric changes in lung volume are represented by a deformed sinusoidal function:

[0079]

[0080] where V r represents functional residual capacity, with units of L, V t represents the baseline tidal volume, with units of L, AV t represents the tidal volume change, R r represents the respiratory rate, with units of bpm.

[0081] During hypoxia, peripheral and central chemoreceptors in the chemical feedback pathway can sense the information of oxygen decrease and carbon dioxide increase. Through neural signals, the central nervous system commands the relevant muscles to increase the respiratory rate and tidal volume, and switch the breathing pattern. The central nervous system model simplifies the transmission of neural signals in the brainstem and designs a detailed piecewise proportional-integral controller, which describes the control of respiratory rate R r and tidal volume change V t through the following equations:

[0082] R r = R r0 + k pr (S r - S e ) + k ir ∫(S r - S e )dt (5)

[0083] AV t = AV t0 + k pv (S r - S e ) + k iv ∫(S r - S e )dt (6)

[0084]

[0085] where S e represents the SpO2 monitored by the chemoreceptors, S r represents the reference value of SpO2, (S r -S e ) represents the deviation value of the feedback signal, k pr represents the proportional controller parameter of respiratory rate, k pv represents the integral control parameter of respiratory rate, k ir represents the proportional control parameter of respiratory amplitude, k iv represents the integral controller parameter of respiratory amplitude, R r0 represents the initial value of respiratory rate, AV t0 represents the initial value of tidal volume change. The variation of controller parameters describes the difference in respiratory control between normal and hypoxic states, and simulates the hyperventilation during hypoxia.

[0086] The cardiovascular system and the respiratory system are two systems that interact closely. The main function of the respiratory system is to maintain blood circulation and material transport in the body through gas exchange using the energy provided by the heart pump. Since the chemoreceptors cannot monitor the oxygen level in real time, they update the information on oxygen demand with each heartbeat. This process is represented as S e = S a (t - τ), where τ represents the delay time of blood circulation.

[0087] In the cardiovascular system model, the heart rate is defined as the number of heartbeats per minute (bpm) and is calculated by 60 / τ, and the heart rate is represented by H r . In terms of blood oxygen and capacity, the oxygen dissociation curve describes the relationship between SpO2 and oxygen partial pressure, which is fitted by the Hill equation, as follows:

[0088]

[0089] where P 50 represents the oxygen partial pressure when the blood oxygen saturation is 50%, with units of mmHg, and r represents the steepness of the oxygen dissociation curve.

[0090] The shape and position of the oxygen dissociation curve are affected by various factors, including body temperature, blood pH, carbon dioxide partial pressure, and the concentration of 2,3-diphosphoglycerate. Under hypoxic conditions, the curve may shift to the left or right.

[0091] The physiological model of human respiration process consists of (1)-(8) and is represented as:

[0092] x = [PA , S a , S v , R r , AV t ] (9)

[0093] u = [P m , H r ] (10)

[0094] y = S a (11)

[0095] Θ = [k1, k2, k3, k4, w, V r , V t , S r , k pr1 , k pr2 ,

[0096] k ir1 , k ir2 , k pv1 , k pv2 , k iv1 , k iv2 , P50, w] (12)

[0097] where x denotes the model state, u denotes the input data, y denotes the output data, and Θ denotes the vector of key parameters.

[0098] From a physiological perspective, the structure of the respiratory model is generally the same for different individuals, while the values of the parameters are individualized. Each experimental subject i has a specific Θ i , called individualized parameters. For the population as a whole, individualized parameters follow a specific distribution. Individualized parameters can be classified according to the adaptability of individuals to hypoxic conditions, which mainly includes two stages: acute hypoxic response and chronic hypoxic response.

[0099] Acute response usually occurs within the first 100 minutes after hypoxic exposure, mainly reflected in changes in heart rate and ventilation. Chronic response occurs within 1-100 days after exposure, mainly reflected in changes in capillary and hemoglobin density. Among them, the increase in capillary density can increase k1 and k2, thereby improving oxygen uptake. Changes in hemoglobin density affect the values of P 50 and r, thereby adjusting the release and utilization efficiency of oxygen. Changes in P 50 reflect the affinity of hemoglobin for oxygen. The higher the r value, the steeper the curve, indicating the stronger the ability of hemoglobin to release oxygen under hypoxic conditions.

[0100] In the chronic hypoxic response, the parameters involved are Θ v = [k1, k2, P 50, r] change over time (in days), these parameters are called day-difference parameters, while other parameters that do not change over days are called fixed parameters, denoted as θ f .

[0101] The second module includes an identifiable-based integrated variational inference (IBC-VI) and hypoxia-aware maximum likelihood estimation (Hypoxia-aware MLE). Establishing the second module includes:

[0102] A joint probability density function is established, which describes the relationship between parameters, input data, and output data, reflecting the hierarchical structure of data. For population and individualized modeling, the joint probability density function can be expressed as:

[0103]

[0104] where θ i represents the individualized parameters of the i-th experimental subject, u i represents the input data of the experimental subject i, y ij represents the j-th measurement value of the experimental subject i;

[0105] The inference target is defined. The inference target is the posterior probability density function, denoted as:

[0106]

[0107] where p(u, y) represents the model evidence, which only depends on the model properties and is independent of the specific parameter estimation results. The higher the model evidence, the better the model fitting degree of the data, and the more reliable the parameter estimation results;

[0108] Variational inference is applied. The principle of variational inference is to find a set of variable parameters v = {v μ , v σ} such that the approximate posterior probability density q(φ, θ|v) is closest to the true posterior probability density p(φ, θ|u, y). The mean vector of the approximate posterior is v μ = [v μ：φ ; v μ：θ ], and the standard deviation vector is v σ = [v σ：φ ; v σ：θ ];

[0109] The optimization target of population and individualized modeling is the evidence lower bound (ELBO). In variational inference, the Kullback-Leibler (KL) divergence is used to measure the difference between the approximate posterior probability density and the true posterior probability density, denoted as:

[0110] DKL (v) = E q [log q(φ, θ | v) - log p(φ, θ | u, y)] (15)

[0111] where D KL (v) represents the KL divergence, E q The [·] operator represents the expectation value. The lower bound of evidence is maximized, denoted as:

[0112] L(v) = log p(y) - D KL (v)

[0113] = E q [log p(y | θ, u) + log p(θ | φ) + log p(φ) - log q(φ, θ | v)] (16)

[0114]

[0115] where L(v) represents the lower bound of evidence, and the lower bound of evidence is calculated by multiple sampling of sample data, and log p(y | θ, u) is represented as:

[0116]

[0117] where l ij represents the length of the signal, represents the simulation output data, represents the real output data. The fitting degree of the model is measured by calculating the difference between the simulation output data and the real output data. The smaller the difference, the larger the likelihood function value, indicating that the model fits the data better.

[0118] log q(φ, θ | v) is represented as:

[0119]

[0120] where ∈ k represents the error of the kth group, v σ represents the precision parameter of the variational distribution of the kth group. The complexity of the approximate posterior distribution is measured by calculating the entropy of the variational distribution, and the larger the entropy, the more complex the approximate posterior distribution.

[0121] The distribution of the individualized parameter θ i is defined, and the individualized parameter θ i obeys a multivariate normal distribution, which is defined by a mean vector and a covariance matrix, represented as:

[0122] θ i ~ N(φ μ , φ ∑) (20)

[0123] where φ μ denotes the mean of the individualized parameters, φ ∑ denotes the correlation between individualized parameters, φ ∑ can be decomposed into a lower triangular matrix φ L and the product of φ L transpose, i.e. Therefore, the variational distribution is represented as φ = {φ μ , φ L};

[0124] The identifiability of parameters is evaluated. Due to the different effects and signal qualities of each experiment, the estimation of parameters θ and φ is not independent. For different experimental objects, if the identifiability of some parameters is better (i.e. more easily and accurately estimated), these high-quality data should be relied on more in the estimation of the parameter distribution of the population model.

[0125] The indicators for evaluating the identifiability of parameters include normalized parameter sensitivity and unnormalized parameter variance. The normalized parameter sensitivity is used to evaluate the identifiability of different parameters in the same experiment, denoted as S i , and the unnormalized parameter variance reflects the identifiability of parameters in experiments of different experimental objects, denoted as V i . The distribution based on the identifiability p(θ|φ) is represented as:

[0126]

[0127] where l θ denotes the dimension of the physiological model parameter vector, and ∑ Vi denotes all first-order and second-order variance matrices of the experimental object i.

[0128] The probability distribution of parameters θ is calculated based on the identifiability of parameters of each experimental object. High-quality data has a higher weight in the estimation of the parameter distribution. The above method allows more accurate estimation of the distribution of model parameters considering the uncertainty and identifiability of parameter estimation.

[0129] An integrated variational inference optimization objective based on identifiability is constructed, which combines the effects of individual model fitting, the rationality of population model, and signal quality, etc. to evaluate the adaptability of the model at the population and individual levels. The probability distributions of parameters θ and φ are represented as log p(θ|φ) and log p(φ), respectively, and the calculation formula is as follows:

[0130]

[0131] log p(φ) = -λ||φ ∑ || *(twenty three)

[0132] Among them, ||.|| * represents the nuclear norm of the matrix, and λ represents a hyperparameter. Substituting Equations (18), (19), (22), and (23) into Equation (16), we obtain the identifiability-based integrated variational inference optimization objective. In this way, the optimization objective can balance the model adaptability at the individual and group levels while taking into account the uncertainty of parameter estimation and data quality, thereby improving the accuracy and robustness of model inference.

[0133] The ensemble variational inference based on identifiability can only calculate the daily average parameters without considering the parameter θ v Changes with the day. Since hypoxemia (SpO2 < 90%) is an important health indicator, the case of SpO2 < 90% should be given a higher weight when estimating the daily difference. In order to take the effect of hypoxemia into account in the optimization process, the hypoxia-aware maximum likelihood estimation (Hypoxia-aware MLE) is used to improve the parameter estimation, which is calculated as follows:

[0134]

[0135] Among them, ∑ w represents a diagonal weight matrix, whose kth element σ w (k) Satisfy:

[0136]

[0137] In the optimization for differentiating daily differences, the parameter thresholds were set to 30% of the full threshold on either side of the mean, so that the estimation of the parameters would take into account a reasonable range of variation around the mean.

[0138] Combining identifiability-based variational inference with hypoxia-aware MLE can simultaneously consider the fitting effect of individual models and the rationality of group models, and capture the parameter θ v The daily changing trend is used to achieve the purpose of parameter estimation.

[0139] The third module includes a sensitivity-based swarm intelligence genetic algorithm (SGBW-GA) and an optimization method based on sensitivity analysis (SA). The third module is used for group and individual parameter estimation. The sensitivity-based swarm intelligence genetic algorithm optimizes the IBC-VI objective function to obtain the θ of each experimental object. f and θ v,ini and φ of the experimental subject population.

[0140] The sensitivity-based swarm intelligence genetic algorithm includes the following steps:

[0141] Step 1, initialize the population, randomly generate multiple individuals containing complete parameter vectors;

[0142] Step 2, according to the order of magnitude of parameter sensitivity, divide the parameters into high sensitivity group g high and low sensitivity group g low ;

[0143] Step 3, calculate the fitness of each individual using the IBC-VI objective function to evaluate its pros and cons;

[0144] Step 4, for high sensitivity group g high , apply genetic algorithm to generate new individuals by selecting individuals with high fitness for crossover and mutation operations;

[0145] Step 5, for low sensitivity group g low , apply genetic algorithm to generate new individuals by selecting individuals with high fitness for crossover and mutation operations;

[0146] Step 6, use set Ω + to select high sensitivity parameters, and extract the gene group of high sensitivity parameters from the chromosome for optimization, and then restore the gene group of high sensitivity parameters to all chromosomes, use set Ω - to select low sensitivity parameters, and extract the gene group of low sensitivity parameters from the chromosome for optimization, and then restore the gene group of low sensitivity parameters to all chromosomes;

[0147] Step 7, repeat the above steps until the preset iteration number is reached or the termination condition is met, to realize the optimization of parameter vector.

[0148] In the sensitivity-based swarm intelligence genetic algorithm, parameter grouping is achieved by selecting the gene group of specific parameters. Use set Ω + to select high sensitivity parameters, and extract the gene group of high sensitivity parameters from the chromosome for optimization, and then restore the gene group of high sensitivity parameters to all chromosomes, use set Ω - to select low sensitivity parameters, and extract the gene group of low sensitivity parameters from the chromosome for optimization, and then restore the gene group of low sensitivity parameters to all chromosomes. The grouping of parameters is based on the order of magnitude of parameter sensitivity ∑ S , where high sensitivity parameters are classified into g high group; low sensitivity parameters are classified into g low group. The sensitivity-based swarm intelligence genetic algorithm realizes the iteration process by alternately optimizing the two groups. Optimizing high sensitivity group helps to accelerate the convergence speed of the algorithm, while optimizing low sensitivity group can record more correct optimization information that is easily overlooked, thereby improving the overall optimization accuracy.

[0149] The sensitivity calculation of the present application adopts Sobol global sensitivity analysis (SA). The sensitivity S describes how to attribute the variance V of the model output to each input variable and its interaction, expressed as:

[0150]

[0151] Wherein, V n represents the first-order variance caused by the input variable θ n , S n represents the corresponding first-order sensitivity, V nm represents the second-order variance caused by the input variable θ n and θ m , S nm represents the corresponding second-order sensitivity, and Var(Y) represents the unconditional variance.

[0152] In the sensitivity-based swarm intelligence genetic algorithm, ∑ S is a diagonal matrix composed of S n . For parameters with the same change amplitude within the threshold value based on fixed input, the higher the parameter sensitivity, the greater the output change. In the identifiable-based integrated variational inference, the matrix ∑ V is composed of V n and V nm .

[0153] Initialize the human respiratory process physiological model in the first module, bring the data in S1 into the second module as real data, use the identifiable-based integrated variational inference and hypoxia perception maximum likelihood estimation for parameter estimation, use the sensitivity-based swarm intelligence genetic algorithm in the third module to optimize the objective function in the second module, and obtain the optimized population and individualized parameters.

[0154] S3, use the test data set to verify the prediction accuracy of the population and individualized human respiratory process physiological model constructed in S2.

[0155] Select a performance evaluation index, and in the present embodiment, the mean square error (MSE) is used to verify the effectiveness of the individualized respiratory model. The calculation method is as follows:

[0156]

[0157] Wherein, l is the signal length of the SpO2 measuring point.

[0158] The performance of the basic genetic algorithm and the sensitivity-based swarm intelligence genetic algorithm in optimizing parameters was compared through experiments. In the case of ensuring the same number of iterations of the two algorithms, multiple experiments were conducted, and the optimal results were compared. A total of 190 cases were analyzed in the experiments, and the SBGW-GA achieved a smaller average MSE in 114 cases, and the basic GA achieved a smaller average MSE in 76 cases. The overall average MSE of the SBGW-GA was 19.19, and the overall average MSE of the basic GA was 19.88. In optimizing parameters, the SBGW-GA generally provides higher accuracy, with an average MSE lower than the basic GA, thus indicating that the SBGW-GA effectively improves the accuracy of the optimization results.

[0159] The simulation results of the population and individualized human respiratory process physiological model of the experimental object 1 are as shown in Figures 1-3 These simulation curves successfully capture the periodic changes in SpO2 during IHT, while reflecting the lung ventilation volume, i.e., the product of the respiratory frequency R r and the tidal volume V t . The SBGW-GA optimization results show that R r increases, while V t slightly decreases, indicating that the central nervous system responds to hypoxia by increasing the respiratory frequency, which is consistent with the actual physiological process. In contrast, the results of the basic GA show that both R r and V t decrease, resulting in a decrease in ventilation volume, which is inconsistent with the physiological phenomenon of excessive ventilation. Therefore, the population and individualized human respiratory process physiological model proposed in the present application can effectively simulate the acute hypoxic response.

[0160] As shown in Figure 4 , two groups of parameter distributions φ were obtained by applying two different optimization objectives, C-VI and IBC-VI. The sampling data of the individualized parameters θ and the parameter distributions φ obtained by the SBGW-GA method, G1 represents the individualized parameters θ, G2 represents the sampling data of φ obtained using the IBC-VI optimization objective, and G3 represents the sampling data of φ obtained using the C-VI optimization objective. Sampling was performed 1000 times, randomly within the theoretical range. The dashed lines in the chart represent the maximum and minimum values of the data, the boundary lines of the box plot represent the first and third quartiles, the median line is located within the box, and the red "+" marks the outliers. The results show that under the two optimization objectives, some parameter distributions are similar, while other parameter distributions differ.

[0161] As shown in Figure 5The population and individualized human respiratory process physiological model simulation curves using different optimization objectives are shown. "Distribution" means 100 random samplings of the population and individualized human respiratory process physiological model. "Mean" means the mean simulation curve of the population and individualized human respiratory process physiological model.

[0162] In the simulation of the two optimization objectives, SpO2 shows periodic oscillation. In the result of IBC-VI optimization objective, the minimum point of hypoxia state gradually rises, which more accurately reflects the gradual adaptation of human body to hypoxia. The simulation curve of IBC-VI optimization objective shows the obvious increase of R r and V t , which more effectively describes the over-ventilation and acute hypoxia reaction than the result of C-VI.

[0163] Overall, the IBC-VI optimization objective helps to correctly identify and exclude low-quality data in the model fitting process, thereby constructing more reasonable population and individualized human respiratory process physiological models.

[0164] The simulation results of the diurnal difference optimization part are shown in Figures 1-3 . Compared with the simulation without considering the diurnal difference parameters, the simulation considering the diurnal difference parameters has higher fitting accuracy, with the average MSE of SpO2 being 15.81.

[0165] Figure 6 The changes of diurnal difference parameters θ v are shown. Each dotted line represents the data curve of an experimental object, and the thick solid line represents the data curves of experimental object 4 and experimental object 6, respectively, which are the experimental objects with the most improved and deteriorated hypoxia tolerance, respectively.

[0166] The shift of oxygen dissociation curve can be further discussed by P 50 and r parameters in Figure 6 . The shift of the curve varies from person to person, which can be left, right or repeated. The change of lung diffusion capacity can be discussed by k1 and k2 parameters in Figure 6 . The greater the increase in lung diffusion capacity, the greater the ΔHTI. The ΔHTI of all experimental objects is distributed in the range of [-78, 309], and the ΔHTI of experimental object 6 is 309, indicating that its hypoxia tolerance is significantly enhanced, and the final values of k1 and k2 are greater than the initial values, which is consistent with the results of the hypoxia tolerance test. The ΔHTI of experimental object 4 is -78, and the k1 and k2 curves show a downward trend, which is also consistent with the experimental results. This shows that by fitting θ v , the population and individualized human respiratory process physiological model can describe the chronic hypoxia reaction.

[0167] Therefore, the application adopts the above-mentioned group and individualized breathing process modeling method based on sensitivity optimization, which realizes accurate simulation of blood oxygen saturation, respiratory frequency and tidal volume, captures hypoxia reaction characteristics, avoids the influence of low-quality signals, and makes the daily difference simulation have higher fitting accuracy.

[0168] Finally, it should be noted that the above examples are only used to illustrate the technical solutions of the present application and not to limit them. Although the present application has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that the technical solutions of the present application can still be modified or equivalently replaced, and these modifications or equivalent replacements cannot make the modified technical solutions deviate from the spirit and scope of the technical solutions of the present application.

Claims

1. A method for modeling the respiratory process of a population and an individual based on sensitivity optimization, characterized in that: The following steps are involved: S1, collect data and preprocess the data; S2. Initialize the group and individualized physiological models of human respiratory process in the first module. The physiological model of human respiratory process is described by the following equation: （9） （10） （11） （12） in, Represents the model state, Represents input data, Indicates output data, represents the key parameter vector, represents the alveolar oxygen partial pressure, Indicates blood oxygen saturation, The diffusion capacity of oxygen from alveoli to arterial blood. Indicates the diffusion capacity of oxygen from arterial blood to tissues. Indicates the diffusion capacity of oxygen from tissues to venous blood. Indicates the diffusion capacity of oxygen from venous blood to alveoli. Indicates venous oxygen saturation, Represents the respiratory rate, Indicates the oxygen partial pressure of the inspired gas, Indicates heart rate, Indicates the effect of metabolism on venous oxygen saturation, Indicates functional residual capacity, represents the basic tidal volume, Indicates the oxygen partial pressure when the blood oxygen saturation is 50%. Indicates the steepness of the oxygen dissociation curve; The data in S1 is brought into the second module as real data, and the parameters are optimized by integrated variational inference based on identifiability and maximum likelihood estimation of hypoxia perception. The optimization objective of integrated variational inference based on identifiability is expressed as: （16） in, Expressed as: （18） Indicates the signal length, represents the simulation output data, Represents the real output data, Represents input data; Expressed as: （19） Indicates the The error of the group, Indicates the The precision parameter of the variational distribution of each component; Expressed as: （22） represents the nuclear norm of the matrix, represents a hyperparameter; Expressed as: （23） represents the mean value of the individualized parameter, represents the correlation between individualized parameters, represents the variational distribution; Add the maximum likelihood estimation of hypoxia perception to the parameter estimation, which is calculated as follows: （24） in, represents the diagonal weight matrix, Indicates the Individualized parameters of each experimental subject, Indicates the experimental subject Input data; The second module is optimized using the sensitivity-based swarm intelligence genetic algorithm in the third module to obtain optimized swarm and individualization parameters; S3. Use the test dataset to verify the prediction accuracy of the group and individualized human respiratory process physiological models constructed in S2.

2. A group and individualized respiratory process modeling system based on sensitivity optimization, applied to the group and individualized respiratory process modeling method based on sensitivity optimization according to claim 1, characterized in that: include: Preprocessing module, used for data collection and preprocessing; The parameter processing module is used to initialize the model and optimize the data to obtain the optimized group and individual parameters.

3. A group and individualized respiratory process modeling system based on sensitivity optimization according to claim 2, characterized in that, The parameter processing module includes: the first module is used to construct group and individualized physiological models of the human respiratory process; the second module uses identifiability-based integrated variational inference and hypoxia perception maximum likelihood estimation for parameter estimation; and the sensitivity-based swarm intelligence genetic algorithm in the third module is used to optimize the objective function in the second module.

4. A group and individual respiratory process modeling system based on sensitivity optimization according to claim 3, characterized in that: The first module includes the respiratory system model, the central nervous system model and the cardiovascular system model.

5. A group and individual respiratory process modeling system based on sensitivity optimization according to claim 4, characterized in that: The central nervous system model simplifies the transmission of neural signals in the brainstem and designs a piecewise proportional-integral controller to describe the control of respiratory rate and tidal volume changes through the following equations: （5） （6） （7） in, Represents the respiratory rate, represents the change in tidal volume, Indicates SpO2 monitored by chemoreceptors, Indicates the reference value of blood oxygen saturation. Indicates the deviation value of the feedback signal, represents the proportional controller parameter of the respiratory frequency, represents the integral control parameter of the respiratory frequency, represents the respiratory amplitude proportional control parameter, represents the integral controller parameter of the breathing amplitude, Indicates the initial value of respiratory rate, Indicates the initial value of tidal volume change.

6. A group and individual respiratory process modeling system based on sensitivity optimization according to claim 5, characterized in that: Daily difference parameters Changes over time.

7. A group and individual respiratory process modeling system based on sensitivity optimization according to claim 3, characterized in that: In the third module, the sensitivity-based swarm intelligence genetic algorithm includes the following steps: Step 1: Initialize the population and randomly generate multiple individuals containing complete parameter vectors; Step 2: According to the order of magnitude of parameter sensitivity, the parameters are divided into high sensitivity groups g high and low sensitivity group g low ; Step 3: Calculate the fitness of each individual using the identifiability-based integrated variational inference objective function to evaluate its quality; Step 4: For high sensitivity group g high Apply genetic algorithms to generate new individuals by selecting individuals with high fitness and performing crossover and mutation operations; Step 5: For low sensitivity group g low Apply genetic algorithms to generate new individuals by selecting individuals with high fitness and performing crossover and mutation operations; Step 6: Use Collections To select high sensitivity parameters, extract the genome of high sensitivity parameters from the chromosome for optimization, and then restore the genome of high sensitivity parameters to all chromosomes, using the set To select low-sensitivity parameters, extract the genome of low-sensitivity parameters from the chromosomes for optimization, and then restore the genome of low-sensitivity parameters to all chromosomes; Step 7: Repeat the above steps until the preset number of iterations is reached or the termination condition is met.

Images