Group and individualized breathing process modeling method based on sensitivity optimization

Abstract

The invention discloses a population and individualized respiratory process modeling method based on sensitivity optimization, and belongs to the technical field of biomedical engineering and physiology, and the method comprises the following steps: S1, data preparation; s2, constructing a model and optimizing parameters; and S3, model verification. According to the population and individualized breathing process modeling method based on sensitivity optimization, accurate simulation of the oxyhemoglobin saturation, the breathing frequency and the tidal volume is achieved, hypoxia reaction characteristics are captured, the influence of low-quality signals is avoided, and daily difference simulation has higher fitting precision.

Description

Technical Field

[0001] The invention belongs to the field of biomedical engineering and physiological technology, and in particular relates to a group and individualized respiratory process modeling method based on sensitivity optimization. Background Art

[0002] Breathing is one of the important physiological functions of the human body. Ensuring normal breathing is the most basic requirement for maintaining life. When people from the plains are in a plateau environment, their breathing will be restricted by the reduced oxygen concentration, which will cause acute altitude sickness and even induce high altitude mountain sickness. Therefore, building a physiological model that can predict indicators such as blood oxygen saturation (SpO2), respiratory rate, and tidal volume is very valuable for the prevention of altitude sickness and hypoxia prep.

[0003] At present, many respiratory models based on physiological principles have been established. The main dynamic characteristics of respiration can be described as the following three aspects. First, pulmonary gas exchange and blood perfusion. Second, the central nervous system (CNS) controls the respiratory process. A CNS model of the respiratory process based on neuronal electrical signal conduction was constructed to accurately describe the respiratory mode and the working mechanism of 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] The construction of group and individual models requires effective parameter estimation schemes. Variational inference (VI) is widely used to evaluate complex damage, and the key lies in selecting a distribution close to the true posterior distribution. Existing technologies optimize the objective through collective variational inference (C-VI) to coordinate low-quality and heterogeneous data, but the signal quality measurement standard has not been established.

[0005] Since existing models face challenges such as limited sample size and inconsistent signal quality, a method to build an individualized physiological model is urgently needed. Summary of the invention

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

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

[0008] S1, collect data and pre-process the data;

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

[0010] S3. Use the test data set to verify the prediction accuracy of the group and individualized human respiratory process physiological models constructed in S2.

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

[0012] Preferably, the parameter processing module includes: a first module for constructing group and individualized physiological models of human breathing process; a second module for performing parameter estimation using integrated variational inference based on identifiability and maximum likelihood estimation of hypoxia perception; and a sensitivity-based swarm intelligence genetic algorithm in a third module for optimizing the objective function in the second module.

[0013] Preferably, the first module includes 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 a piecewise proportional-integral controller is designed to describe the control of respiratory rate 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] Among them, R r represents the respiratory rate, ΔVt Indicates the change in tidal volume, S e represents SpO2 monitored by chemoreceptors, S r Indicates the reference value of blood oxygen saturation, (S r -S e ) represents the deviation value of the feedback signal, k pr represents the proportional controller parameter of the respiratory frequency, k pv represents the integral control parameter of the respiratory frequency, k ir represents the respiratory amplitude proportional control parameter, k iv represents the integral controller parameter of the breathing amplitude, R r0 Indicates the initial value of respiratory rate, ΔV t0 Indicates the initial value of tidal volume change.

[0019] Preferably, the physiological model of human breathing process is described by the following equation:

[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] Among them, x represents the model state, u represents the input data, y represents the output data, θ represents the key parameter vector, and P A represents alveolar oxygen partial pressure, S a represents blood oxygen saturation, k1 represents the diffusion capacity of oxygen from alveoli to arterial blood, k2 represents the diffusion capacity of oxygen from arterial blood to tissues, k3 represents the diffusion capacity of oxygen from tissues to venous blood, k4 represents the diffusion capacity of oxygen from venous blood to alveoli, S vrepresents venous oxygen saturation, R r Respiratory rate, P m Indicates the oxygen partial pressure of the inspired gas, H r represents heart rate, w represents the effect of metabolism on venous oxygen saturation, V r Indicates functional residual capacity, V t represents the basic tidal volume, P 50 It represents the oxygen partial pressure when the blood oxygen saturation is 50%, and r represents the steepness of the oxygen dissociation curve.

[0026] Preferably, the daily difference parameter θ v =[k1, k2, P 50 , r] changes over time.

[0027] Preferably, in the second module, the optimization objective of the integrated variational reasoning based on identifiability is expressed 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] Where log p(y|θ,u) is expressed as:

[0031]

[0032] l ij Indicates the signal length, represents the simulation output data, represents the real output data, and u represents the input data;

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

[0034]

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

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

[0037]

[0038] ||.|| * represents the nuclear norm of the matrix, and λ represents the hyperparameter;

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

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

[0041] φ μ represents the average value of the individualized parameter, φ Σ represents the correlation between individualized parameters, and φ represents the variational distribution.

[0042] Preferably, in the second module, the maximum likelihood estimation of hypoxia perception is added to the parameter estimation, and the calculation method is as follows:

[0043]

[0044] Among them, ∑ w represents the diagonal weight matrix, θ i represents the individualized parameter of the ith experimental subject, u i Represents the input data of experimental subject i.

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

[0046] Step 1: Initialize the population and randomly generate multiple individuals containing complete parameter vectors;

[0047] 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 ;

[0048] Step 3: Calculate the fitness of each individual using the identifiability-based integrated variational inference objective function to evaluate its quality;

[0049] 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;

[0050] Step 5: Low sensitivity group g low Apply genetic algorithm to 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, extract the genome of high-sensitivity parameters from the chromosomes 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;

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

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

[0054] (1) The present invention incorporates the control of the nervous system over the respiratory system into the modeling, thereby achieving the simulation of indicators such as blood oxygen saturation, respiratory rate, and tidal volume, and capturing the characteristics of acute and chronic hypoxic responses;

[0055] (2) The design of the integrated variational inference (IBC-VI) optimization objective based on identifiability in this invention innovatively measures the signal quality, avoids the negative impact of low-quality signals on distribution estimation, and improves the adaptability of the multi-level evaluation model for groups and individuals;

[0056] (3) The present invention proposes an optimization scheme based on sensitivity analysis (SA) and a swarm intelligence genetic algorithm based on sensitivity (SGBW-GA), which effectively improves the optimization accuracy by calculating the sensitivity of model parameters and optimizing the parameters with large differences respectively;

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

[0058] (5) The daily difference simulation of the present invention has higher fitting accuracy, which can more clearly understand the complex physiological processes of the human body and analyze the similarities and differences between people.

[0059] The technical solution of the present invention is further described in detail below through the accompanying drawings and embodiments. BRIEF DESCRIPTION OF THE DRAWINGS

[0060] Figure 1 This is a SpO2 data graph of a group and individualized respiratory process modeling method based on sensitivity optimization of the present invention. 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 parameter for distinguishing daily differences. 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 parameter for distinguishing daily differences 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 2It is a respiratory frequency curve diagram of a group and individualized respiratory process modeling method based on sensitivity optimization of the present invention. Each sub-diagram includes three groups. The three groups represent the real curve, the SBGW-GA simulation curve and the basic GA simulation curve from top to bottom, respectively. a represents the change of respiratory frequency of the three different groups over time on the 2nd day; b represents the change of respiratory frequency of the three different groups over time on the 6th day; c represents the change of respiratory frequency of the three different groups over time on the 10th day;

[0062] Figure 3 It is a tidal volume curve diagram of a group and individualized respiratory process modeling method based on sensitivity optimization of the present invention. Each sub-diagram includes three groups. The three groups respectively represent the basic GA simulation curve, the curve of daily difference parameter result simulation and the SBGW-GA simulation curve from top to bottom. a represents the change of tidal volume of the three different groups over time on the second day; b represents the change of tidal volume of the three different groups over time on the sixth day; c represents the change of tidal volume of the three different groups over time on the tenth day;

[0063] Figure 4 The box plot of individualized parameters and distribution sampling of a group and individualized respiratory process modeling method based on sensitivity optimization of the present invention includes three different groups in each sub-plot, and these three groups are individualized parameters θ from left to right. i , the sampling data of φ obtained using the IBC-VI optimization target and the sampling data of φ obtained using the C-VI optimization target, a represents the difference in parameter k1 between different groups; b represents the difference in parameter k2 between different groups; c represents the difference in parameter P between different groups 50 d represents the difference of parameter r between different groups, and “+” represents an outlier;

[0064] Figure 5 It is a physiological parameter simulation curve diagram of a group and individualized respiratory process modeling method based on sensitivity optimization of the present invention, wherein a shows the change of blood oxygen saturation at different time points, and the mean value of IBC-VI gradually increases; b shows the change of tidal volume over time, and as time goes by, the mean value of tidal volume gradually increases, and the mean value of IBC-VI is always lower than C-VI; c shows the change of respiratory frequency over time, and as time goes by, the mean value of respiratory frequency also gradually increases, and the mean value of IBC-VI is always higher than C-VI;

[0065] Figure 6 is the daily difference parameter θ of a group and individual respiratory process modeling method based on sensitivity optimization of the present invention. vThe graph 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 solution and advantages of the embodiments of the present invention clearer, the technical solution 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 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 meaning understood by people with general skills in the field to which the present invention belongs.

[0067] Embodiment 1

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

[0069] S1, using 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 oxygen partial pressure of 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 periodic 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 days of IHT, reflecting the changes in the subjects' hypoxia adaptability. The changes in hypoxia tolerance indicators are 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: a first module for constructing group and individualized physiological models of human respiratory process; a second module for optimizing parameters using integrated variational reasoning based on identifiability and maximum likelihood estimation of hypoxia perception; and a third module for optimizing the second module using a sensitivity-based group intelligence genetic algorithm.

[0073] The first module is used to construct a group and individualized physiological model 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 arterial blood vessels. However, in a hypoxic environment, the diffusion of gases 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 breathing ventilation in the first few minutes when a person is exposed to a hypoxic environment. The response will also weaken after a few tens of minutes after exposure. 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 It represents the 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 vThe derivative, w represents the effect of metabolism on SvO2, the unit is 1 / s, k1 represents the diffusion capacity of oxygen from alveoli to arterial blood, k2 represents the diffusion capacity of oxygen from arterial blood to tissues, k3 represents the diffusion capacity of oxygen from tissues to venous blood, k4 represents the diffusion capacity of oxygen from venous blood to alveoli, k1 unit is mmHg / s, k2, k3 and k4 unit is 1 / s, f d Represents the oxygen dissociation curve function.

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

[0079]

[0080] Among them, V r Indicates functional residual gas volume, in L, V t Indicates the basic tidal volume, in L, ΔV t Indicates the change in tidal volume, R r Indicates the respiratory rate in bpm.

[0081] During hypoxia, peripheral and central chemical receptors in the chemical feedback pathway sense the drop in oxygen and the rise in carbon dioxide. Through neural signals, the central nervous system commands the relevant muscles to increase the respiratory rate and tidal volume and switch the breathing mode. The central nervous system model simplifies the transmission of neural signals in the brainstem and designs a detailed piecewise proportional-integral controller to describe the respiratory rate R through the following equation r and tidal volume changes V t Control:

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

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

[0084]

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

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

[0087] In the cardiovascular system model, heart rate is defined as the number of beats per minute (bpm) and is calculated by 60 / τ. r In terms of blood oxygenation 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] Among them, P 50 It indicates the oxygen partial pressure when the blood oxygen saturation is 50%, in mmHg, and r indicates the steepness of the oxygen dissociation curve.

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

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

[0092] x=[PA , S a , S v , R r , ΔV 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] Among them, x represents the model state, u represents the input data, y represents the output data, and θ represents the key parameter vector.

[0098] From a physiological point of view, the structure of the breathing model is usually the same for different individuals, but the values ​​of the parameters are individual. Each experimental subject i has a specific θ i , called individualized parameters. For the whole, individualized parameters follow a specific distribution. Individualized parameters can be classified according to the individual's adaptive changes to hypoxic conditions. The adaptive changes are mainly divided into two stages: acute hypoxic response and chronic hypoxic response.

[0099] Acute reactions usually occur within the first 100 minutes after hypoxia exposure, mainly manifested by changes in heart rate and ventilation. Chronic reactions occur within 1-100 days after exposure, mainly manifested by changes in capillary and hemoglobin density. Among them, the increase in capillary density may increase k1 and k2, thereby improving oxygen uptake. Changes in hemoglobin density will affect P 50 and r, thereby adjusting the release and utilization efficiency of oxygen. 50 The change in reflects the affinity of hemoglobin for binding to oxygen. The higher the r value, the steeper the curve, indicating that hemoglobin has a stronger ability to release oxygen under hypoxic conditions.

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

[0101] The second module includes integrated variational inference based on identifiability (IBC-VI) and hypoxia-aware maximum likelihood estimation (Hypoxia-aware MLE). The second module is built by:

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

[0103]

[0104] Among them, θ i represents the individualized parameter of the ith experimental subject, u i represents the input data of experimental subject i, y ij represents the jth measurement value of experimental object i;

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

[0106]

[0107] Among them, p(u, y) represents model evidence, which depends only on the model properties and not on the specific parameter estimation results. The higher the model evidence, the better the model fits the data and the more reliable the parameter estimation results;

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

[0109] The optimization goal of group and individual 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, expressed as:

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

[0111] Among them, D KL (v) represents KL divergence, E q The [·] operator indicates the expected value. The maximization operation of the lower limit of evidence is expressed 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] Among them, L(v) represents the lower limit of evidence, which is calculated by sampling data for multiple times. Log p(y|θ,u) is expressed as:

[0116]

[0117] Among them, l ij Indicates the signal length, represents the simulation output data, Represents the real output data. The difference between the simulated output data and the real output data is calculated to measure the degree of fit of the model. The smaller the difference, the larger the likelihood function value, indicating that the model fits the data better.

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

[0119]

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

[0121] Define the individualization parameter θ i The distribution of individual parameters θ i It obeys the multivariate normal distribution, which is defined by the mean vector and the covariance matrix, expressed as:

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

[0123] Among them, φ μ represents the average value of the individualized parameter, φ ∑ represents the correlation between individualized parameters, φ ∑ It can be decomposed into the lower triangular matrix φ L and φ L The product of the transpose, that is Therefore, the variational distribution is expressed as φ = {φ μ ,φ L};

[0124] Assess the identifiability of parameters. Since the effects and signal quality of each experiment are different, the estimates of the parameters θ and φ are not independent. If some parameters are more identifiable (i.e., easier to estimate accurately) for different experimental subjects, then more reliance should be placed on these high-quality data in the estimation of the parameter distribution of the population model.

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

[0126]

[0127] Among them, l θ represents the dimension of the physiological model parameter vector, ∑ Vi represents all first-order and second-order variance matrices of experimental subject i.

[0128] The probability distribution of the parameters θ is calculated based on the identifiability of the parameters for each experimental subject. High-quality data has a higher weight in the estimation of parameter distributions. The above method allows for a more accurate estimation of the distribution of model parameters while taking into account the uncertainty and identifiability of parameter estimates.

[0129] An identifiability-based integrated variational inference optimization objective is constructed, which combines factors such as individual model fitting effect, group model rationality, and signal quality to evaluate the adaptability of the model at the group and individual levels. The probability distribution of parameters θ and φ is represented by 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 the hyperparameter. Substituting formula (18), formula (19), formula (22), and formula (23) into formula (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, a higher weight should be given to the case of SpO2 < 90% when estimating the daily difference. In order to consider the effect of hypoxemia 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 distinguishing 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 and randomly generate multiple individuals containing complete parameter vectors;

[0142] 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 ;

[0143] Step 3: Use the IBC-VI objective function to calculate the fitness of each individual to evaluate its quality;

[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: 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 the collection... + To select high-sensitivity parameters, extract the genome of high-sensitivity parameters from the chromosomes for optimization, and then restore the genome of high-sensitivity parameters to all chromosomes, using the collection… - 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;

[0147] Step 7: Repeat the above steps until the preset number of iterations is reached or the termination condition is met to achieve the optimization of the parameter vector.

[0148] In the sensitivity-based swarm intelligence genetic algorithm, parameter grouping is achieved by selecting a genome of specific parameters. + To select high-sensitivity parameters, extract the genome of high-sensitivity parameters from the chromosomes 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 chromosome for optimization, and then restore the genome of low-sensitivity parameters to all chromosomes. The grouping basis of parameters is parameter sensitivity ∑ S of the order of magnitude, where the high sensitivity parameter is included in g high group; parameters with low sensitivity are grouped into g low The sensitivity-based swarm intelligence genetic algorithm implements the iterative process by alternately optimizing the two groups. Optimizing the high-sensitivity group helps accelerate the convergence of the algorithm, while optimizing the low-sensitivity group can record more correct optimization information that is easily overlooked, thereby improving the overall optimization accuracy.

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

[0150]

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

[0152] In the swarm intelligence genetic algorithm based on sensitivity, ∑ S By S n For fixed input and parameters with the same change amplitude within the threshold, the higher the parameter sensitivity, the greater the output change. In identifiability-based integrated variational inference, the matrix ∑ V By V n and V nm composition.

[0153] Initialize the physiological model of human breathing process in the first module, bring the data in S1 into the second module as real data, use integrated variational reasoning based on identifiability and maximum likelihood estimation of hypoxia perception for parameter estimation, and use the swarm intelligence genetic algorithm based on sensitivity in the third module to optimize the objective function in the second module to obtain the optimized swarm and individualized parameters;

[0154] S3. Use the test data set to verify the prediction accuracy of the group and individualized human respiratory process physiological models constructed in S2.

[0155] The performance evaluation index is selected. In this embodiment, the mean square error (MSE) is used to verify the effectiveness of the individualized breathing model. The calculation method is as follows:

[0156]

[0157] Where l is the signal length of the SpO2 measurement point.

[0158] The performance of the basic genetic algorithm and the sensitivity-based swarm intelligence genetic algorithm in optimizing parameters were experimentally compared. Multiple experiments were conducted while ensuring that the number of iterations of the two algorithms was the same, and the best results were compared. A total of 190 cases were analyzed in the experiment, and the SBGW-GA achieved a smaller average MSE in 114 cases, while 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. The SBGW-GA provided higher accuracy overall when optimizing parameters, and its average MSE was lower than that of the basic GA, thus indicating that the SBGW-GA effectively improved the accuracy of the optimization results.

[0159] The simulation results of the physiological model of the group and individual human respiratory process of experimental subject 1 are as follows Figure 1-Figure 3 These simulation curves successfully capture the periodic changes of SpO2 during IHT and reflect the pulmonary ventilation, that is, the respiratory rate R r and tidal volume V t The SBGW-GA optimization results show that R r Increase, while V t The results of basic GA showed that R r and V t All of them are decreasing, resulting in a decrease in ventilation volume, which is inconsistent with the physiological phenomenon of hyperventilation. Therefore, the group and individualized human respiratory process physiological model proposed by the present invention can effectively simulate the acute hypoxic response.

[0160] like Figure 4 As shown in the figure, two sets of parameter distribution φ are obtained by applying two different optimization objectives, C-VI and IBC-VI. The sampling data of individualized parameters θ and parameter distribution φ obtained by the SBGW-GA method, G1 represents the individualized parameter θ, 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. The sampling was performed 1000 times, randomly within the theoretical range. The dotted 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, the distributions of some parameters are similar, while the distributions of other parameters are different.

[0161] like Figure 5As shown in the figure, the group and individual human respiratory process physiological model simulation curves of parameter distribution φ obtained using different optimization objectives. "Distribution" means 100 random samplings of the group and individual human respiratory process physiological models. "Mean" represents the mean simulation curves of the group and individual human respiratory process physiological models.

[0162] In the simulation of the two optimization targets, SpO2 showed periodic oscillations. In the results of the IBC-VI optimization target, the lowest point of the hypoxic state gradually increased, which more accurately reflects the gradual adaptation process of the human body to hypoxia. The simulation curve of the IBC-VI optimization target shows that R r and V t This is a clear increase in the number of systolic and diastolic systolic pressures, which more effectively describes the hyperventilation and acute hypoxic response than the C-VI results.

[0163] In general, the IBC-VI optimization target helps to correctly identify and exclude low-quality data during the model fitting process, thereby constructing a more reasonable group and individualized physiological model of the human respiratory process.

[0164] The simulation results of the daily difference optimization part are as follows: Figure 1-Figure 3 As shown in Figure 2, the simulation with daily difference parameters taken into account has higher fitting accuracy than the simulation without considering daily difference parameters, where the average MSE of SpO2 is 15.81.

[0165] Figure 6 The daily difference parameter θ is shown v , each dotted line represents the data curve of an experimental subject, and the thick solid lines represent the data curves of experimental subject 4 and experimental subject 6, which are the experimental subjects with the most improved and deteriorated hypoxia tolerance, respectively.

[0166] The shift in the oxygen dissociation curve can be determined by Figure 6 P 50 The curve shift varies from person to person and may shift to the left, right, or reciprocally. Changes in lung diffusion capacity can be measured by hypoxia tolerance tests and Figure 6 The more the lung diffusion capacity is enhanced, the greater the ΔHTI. The ΔHTI of all experimental subjects is distributed in the range of [-78,309]. The ΔHTI of experimental subject 6 is 309, indicating that its hypoxia tolerance is significantly enhanced. 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 subject 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 , group and individualized physiological models of human respiratory processes can describe chronic hypoxic responses.

[0167] Therefore, the present invention adopts the above-mentioned group and individualized respiratory process modeling method based on sensitivity optimization, which realizes accurate simulation of blood oxygen saturation, respiratory rate and tidal volume, captures the characteristics of hypoxia response, 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 embodiments are only used to illustrate the technical solution of the present invention rather than to limit it. Although the present invention has been described in detail with reference to the preferred embodiments, those skilled in the art should understand that they can still modify or replace the technical solution of the present invention with equivalents, and these modifications or equivalent replacements cannot cause the modified technical solution to deviate from the spirit and scope of the technical solution of the present invention.

Claims

1. A method for modeling a population and individualized respiratory process based on sensitivity optimization, characterized in that: The following steps are involved: S1, collect data and pre-process the data; S2, initialize the group and individualized physiological models of human respiratory process in the first module, bring the data in S1 into the second module as real data, optimize the parameters based on integrated variational reasoning of identifiability and maximum likelihood estimation of hypoxia perception, and use the sensitivity-based group intelligence genetic algorithm in the third module to optimize the second module to obtain the optimized group and individualized parameters; S3. Use the test data set 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, 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: a first module for constructing group and individualized physiological models of human respiratory process; a second module for parameter estimation using integrated variational inference based on identifiability and maximum likelihood estimation of hypoxia perception; and a sensitivity-based swarm intelligence genetic algorithm in the third module for optimizing the objective function in the second module.

4. A group and individualized 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 individualized 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: R r =R r0 +k pr (S r -S e )+k ir ∫(S r -S e )dt (5) ΔV t =ΔV t0 +k pv (s r -S e )+k iv ∫(s r -S e )dt (6) Among them, R r represents the respiratory rate, ΔV t Indicates the change in tidal volume, S e represents SpO2 monitored by chemoreceptors, S r Indicates the reference value of blood oxygen saturation, (S r -S e ) represents the deviation value of the feedback signal, k pr represents the proportional controller parameter of the respiratory frequency, k pv represents the integral control parameter of the respiratory frequency, k ir represents the respiratory amplitude proportional control parameter, k iv represents the integral controller parameter of the breathing amplitude, R r0 Indicates the initial value of respiratory rate, ΔV t0 Indicates the initial value of tidal volume change.

6. A group and individualized respiratory process modeling system based on sensitivity optimization according to claim 3, characterized in that: The physiological model of human breathing process is described by the following equation: x=[P A ,S a ,S v ,R r ,ΔV t ] (9) u=[P m ,H r ] (10) and=S a (11) θ=[k1,k2,k3,k4,w,V r ,V t ,S r ,k pr1 ,k pr2 ,k ir1 ,k ir2 ,k pv1 ,k pv2 ,k iv1 ,k iv2 ,P 50 ,r] (12) Among them, x represents the model state, u represents the input data, y represents the output data, θ represents the key parameter vector, and P A represents alveolar oxygen partial pressure, S a represents blood oxygen saturation, k1 represents the diffusion capacity of oxygen from alveoli to arterial blood, k2 represents the diffusion capacity of oxygen from arterial blood to tissues, k3 represents the diffusion capacity of oxygen from tissues to venous blood, k4 represents the diffusion capacity of oxygen from venous blood to alveoli, S v represents venous oxygen saturation, R r Respiratory rate, P m Indicates the oxygen partial pressure of the inspired gas, H r represents heart rate, w represents the effect of metabolism on venous oxygen saturation, V r Indicates functional residual capacity, V t represents the basic tidal volume, P 50 It represents the oxygen partial pressure when the blood oxygen saturation is 50%, and r represents the steepness of the oxygen dissociation curve.

7. A group and individualized respiratory process modeling system based on sensitivity optimization according to claim 6, characterized in that: Daily difference parameter θ v =[k1,k2,P 50 ,r] changes over time.

8. A group and individualized respiratory process modeling system based on sensitivity optimization according to claim 3, characterized in that: In the second module, the optimization objective of identifiability-based integrated variational inference is expressed as: L(v)=log p(y)-D KL (v) =E q [log p(y|θ,u)+log p(θ|φ)+log p(φ)-log q(φ,θ|v)] (16) Where log p(y|θ,u) is expressed as: l ij Indicates the signal length, represents the simulation output data, represents the real output data, and u represents the input data; log q(φ,θ|v) is expressed as: ∈ k represents the error of the kth group, v σ represents the precision parameter of the variational distribution of the kth component; log p(θ|φ) is expressed as: ||.|| * represents the nuclear norm of the matrix, and λ represents the hyperparameter; log p(φ) is expressed as: log p(φ)=-λ||φ ∑ || * (23) φ μ represents the average value of the individualized parameter, φ Σ represents the correlation between individualized parameters, and φ represents the variational distribution.

9. A group and individualized respiratory process modeling system based on sensitivity optimization according to claim 3, characterized in that: In the second module, the maximum likelihood estimation of hypoxia perception is added to the parameter estimation, and the calculation method is as follows: Among them, ∑ w represents the diagonal weight matrix, θ i represents the individualized parameter of the ith experimental subject, u i Represents the input data of experimental subject i.

10. A group and individualized 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 algorithm to generate new individuals by selecting individuals with high fitness for crossover and mutation operations; Step 5: Low sensitivity group g low Apply genetic algorithm to generate new individuals by selecting individuals with high fitness for crossover and mutation operations; Step 6: Use the set Ω + To select high-sensitivity parameters, extract the genome of high-sensitivity parameters from the chromosomes 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