Method and device for determining drug combination and storage medium

Abstract

The invention discloses a method for determining drug combination, which comprises the following steps: constructing a calculation model, and constructing the calculation model according to a response surface method; obtaining a reaction value, and inputting the reaction value to the calculation model, the reaction value including a factor main effect and a multi-factor interaction effect; a regression algorithm adapted to the relationship between the reaction variable and the reaction value; the optimal solution of the drug combination is obtained, the obtaining mode of the optimal solution is obtained from the regression algorithm, and the method further comprises a corresponding device and a storage medium. The method has the advantages that the calculation model for drug combination is constructed through the response surface method, the optimal solution of drug combination is obtained through calculation by inputting the factors, and compared with an exhaustion mode in the prior art, the efficiency is higher, and more time, manpower and resource consumption are saved.

Description

Technical Field

[0001] The present invention relates to the field of assisted medical technologies, and particularly to a method, apparatus, and storage medium for determining combined medication. Background Art

[0002] In the treatment of tumor diseases, combined medication is widely adopted as a medication method. Specifically, multiple different drugs are combined and applied, and each drug has its own systemic toxicity, treatment mechanism, etc. It should be noted that each drug can better exert its efficacy only when combined in a specific proportion. Therefore, it involves dosage or dosing factors to determine the optimal dispensing result of the drug.

[0003] Based on the general strategy in mathematics, for the problem of determining the dosage of each of the above drugs, by means of exhaustive enumeration, the operation is very cumbersome. For example, in the case of 4 drugs, each having 3 dosage levels, by means of exhaustive enumeration, 81 operations are required. It can be seen that the time, manpower, and resource consumption are serious, and the efficiency is very low. Summary of the Invention

[0004] In view of the above-mentioned defects, the purpose of the present invention is to provide a method for determining combined medication.

[0005] The technical solution of the present invention is as follows:

[0006] A method for determining combined medication, comprising:

[0007] Constructing a calculation model, and the calculation model is constructed according to the response surface method;

[0008] Obtaining a response value, and inputting the response value into the calculation model, the response value including factor main effects and multi-factor interaction effects;

[0009] Fitting a linear regression algorithm and / or a non-linear regression algorithm suitable for the relationship between the response variable and the response value;

[0010] Obtaining the optimal solution of combined medication, and the optimal solution is obtained by means of the linear regression algorithm and / or the non-linear regression algorithm.

[0011] In some preferred embodiments, the step of fitting a linear regression algorithm and / or a non-linear regression algorithm suitable for the relationship between the response variable and the response value includes:

[0012] Fitting a first-order regression model and / or a second-order regression model;

[0013] Inputting the response value into the first-order regression model and / or the second-order regression model, and obtaining variable values;

[0014] An algorithm for checking the inadequacy of fitting the first-order regression model or / and the second-order regression model;

[0015] Input the variable values into the inadequacy checking algorithm and obtain the checking result;

[0016] Input the checking result and make an adaptability judgment on the first-order regression model or / and the second-order regression model according to the checking result;

[0017] Obtain the first-order regression model or / and the second-order regression model that meets the adaptability.

[0018] In some preferred embodiments, it includes:

[0019] Calculate the stationary point value, and the stationary point value is calculated by the first-order regression model or / and the second-order regression model that meets the adaptability;

[0020] Judge whether the stationary point value meets the optimal solution for combined medication. If so, adopt the corresponding stationary point value.

[0021] In some preferred embodiments, the step of obtaining the optimal solution for combined medication, and the obtaining method of the optimal solution includes obtaining from the linear regression algorithm or / and the non-linear regression algorithm and further includes:

[0022] Construct a hybrid optimization strategy algorithm, and calculate the calculated value for obtaining the optimal solution for combined medication through this hybrid optimization strategy algorithm to obtain the maximum efficacy value;

[0023] Set the constraint conditions for the calculation method, and the constraint conditions include one or more of the single-dose range, the daily-dose range, and the safety limit.

[0024] In some preferred embodiments, the hybrid optimization strategy algorithm includes one of the genetic algorithm, the simulated annealing algorithm, and the particle swarm algorithm, and the simulated annealing algorithm is used for global search, and the particle swarm algorithm is used for local optimization.

[0025] In some preferred embodiments, the step of constructing the calculation model, and the calculation model is constructed according to the response surface method and includes:

[0026] Construct a single-drug PK calculation equation;

[0027] Construct a PK correction calculation equation for the single-drug PK calculation equation;

[0028] Calculate the AUC value based on the above single-drug PK calculation equation,

[0029] Construct a response surface calculation model based on the above PK correction and AUC value.

[0030] In some preferred embodiments, the step of constructing the calculation model, and the calculation model constructed according to the response surface method includes:

[0031] Determine the numerical value of the drug dosage level;

[0032] Input the patient's individualized parameter data;

[0033] Establish an orthogonal test table.

[0034] In some preferred embodiments, the individualized parameter data includes: basic physiological parameters, genetic factors, and disease status.

[0035] A device for implementing the method for determining combined drug use according to any one of the above, and includes: a model architecture module for constructing a calculation model, an experimental design module for constructing an experimental plan, an optimization solution module for calculating the optimal solution of combined drug use, and a distributed calculation module for accelerating the calculation process.

[0036] A readable storage medium, on which a program or instruction is stored, and when the program or instruction is executed by a processor, the steps of the method according to any one of the above are implemented.

[0037] Advantageous effects: The present invention constructs a calculation model for combined drug use by the response surface method, and by inputting factors and calculating the optimal solution of combined drug use. Compared with the exhaustive method in the prior art, it has higher efficiency, saves more time, manpower, and resource consumption; in addition, the compartment model equation in pharmacokinetics is introduced in the response surface calculation to realize accurate modeling of dynamic drug effects, significantly improving the scientificity of combined drug use and achieving precise optimization. Compared with simply superimposing the effects of two drugs in the prior art, only evaluating the drug effect based on the administered dose, the innovative modeling method based on AUC overcomes the limitations of simple dose indicators, provides a more comprehensive assessment of the drug exposure-effect relationship, takes into account factors such as age, liver and kidney function, and gene polymorphism, significantly reduces the bleeding risk, optimizes the test efficiency, and reduces the sample collection volume. BRIEF DESCRIPTION OF THE DRAWINGS

[0038] Figure 1 It is a schematic diagram of the main steps of the method for determining combined drug use in the present invention.

[0039] Figure 2 It is a schematic diagram of a sub-step in the fitting calculation of the method for determining combined drug use in the present invention.

[0040] Figure 3 It is a schematic diagram of another sub-step in the fitting calculation of the method for determining combined drug use in the present invention.

[0041] Figure 4 It is a schematic diagram of another sub-step in the fitting calculation of the method for determining combined drug use in the present invention.

[0042] Figure 5 It is a schematic diagram of another sub-step in the fitting calculation of the method for determining the combined use of drugs in the present invention.

[0043] Figure 6 It is a schematic diagram of a sub-step in constructing a calculation model for the method for determining the combined use of drugs in the present invention.

[0044] Figure 7 It is a schematic diagram of another sub-step in constructing a calculation model for the method for determining the combined use of drugs in the present invention.

[0045] Figure 8 It is a schematic diagram of the structure of the device in the present invention. Detailed implementation manners

[0046] The present invention will be further described below in conjunction with the accompanying drawings and embodiments:

[0047] As Figure 1 shown, a method for determining the combined use of drugs includes the following steps:

[0048] S100. Construct a calculation model, and the calculation model is constructed according to the response surface method;

[0049] S200. Obtain response values and input the response values into the calculation model, where the response values include factor main effects and multi-factor interaction effects;

[0050] S300. Fit a linear regression algorithm and / or a non-linear regression algorithm suitable for the relationship between the response variable and the response values;

[0051] S400. Obtain the optimal solution for the combined use of drugs, and the optimal solution is obtained by means of obtaining from the linear regression algorithm and / or the non-linear regression algorithm.

[0052] Furthermore, as Figure 2 shown, the regression algorithm for fitting the relationship between the response variable and the response values in the above steps includes:

[0053] S301. Fit a first-order regression model;

[0054] S302. Input the response values into the first-order regression model;

[0055] S303. Fit the lack-of-fit test algorithm for the first-order regression model;

[0056] S304. Input the variable values into the lack-of-fit test algorithm and obtain the test result;

[0057] S305. Input the test result and make an adaptability judgment on the first-order regression model according to the test result;

[0058] S306. Obtain a first-order regression model that meets the adaptability.

[0059] In step S301, the fitting of the first-order regression model can be established by the least squares method for pattern fitting to establish a first-order fitting equation for calculation.

[0060] In step S302, the response value is measured through actual clinical examinations. For example, if the goal is to reduce the patient's liver function indicators, the response value can be the AST index value in liver function, and the input data can be the value after administration of a drug combination (drugs a, b, and c).

[0061] Steps S303 - S306 aim to ultimately obtain a first-order regression model that meets the requirements. If a first-order regression model that meets the adaptability cannot be obtained in the judgment of step S305, steps S301 - 305 are repeated.

[0062] When fitting the first-order regression model, the method of steepest ascent and the method of steepest descent are introduced for iterative calculation until the first-order regression model is no longer suitable and there is curvature in the system, and then the second-order regression model is further adapted.

[0063] Furthermore, as Figure 3 shown, the steps are applicable to the regression algorithm for the relationship between the response variable and the response value, including:

[0064] S310. Fit a second-order regression model;

[0065] S320. Input the response value into the second-order regression model;

[0066] S330. Fit the lack-of-fit test algorithm for the second-order regression model;

[0067] S340. Input the variable value into the lack-of-fit test algorithm and obtain the test result;

[0068] S350. Input the test result and make an adaptability judgment on the second-order regression model according to the test result;

[0069] S360. Obtain a second-order regression model that meets the adaptability.

[0070] In step S310, the fitting of the second-order regression model is established by the fitting experiment of the Central Composite Design (CCD) for the experiment design of fitting the second-order regression model.

[0071] In step S320, the response value is measured through actual clinical examinations. For example, if the goal is to reduce the patient's liver function indicators, the response value can be the numerical value of the AST indicator in liver function, and the input data can be the numerical values after administering a drug combination (drugs a, b, and c).

[0072] Steps S330 - S360 aim to finally obtain a first-order regression model that meets the requirements. If a second-order regression model that meets the adaptability cannot be obtained in the judgment of step S305, steps S301 - 305 are repeated.

[0073] In the use of the second-order regression model, the best solution for combined drug use is obtained from the calculated numerical stationary point. Determine the location of this numerical stationary point as a maximum, minimum, or saddle point. If the location of the numerical stationary point is a maximum or minimum, it is determined as the best solution. If the location of the numerical stationary point is a saddle point, it is determined as a value with room for improvement, and then the above steps are repeated.

[0074] Comprehensively, as Figure 4 shown, at the end of the processing through the first-order and second-order regression models, the following steps are included:

[0075] S300 - 1: Calculate the numerical value of the stationary point, and the numerical value of the stationary point is calculated through the first-order regression model or / and second-order regression model that meets the adaptability;

[0076] S300 - 2: Determine whether the numerical value of the stationary point meets the best solution for combined drug use. If so, adopt the corresponding numerical value of the stationary point.

[0077] In specific implementation, the above first-order regression model fitting and second-order regression model fitting can be carried out sequentially, that is, first fit out the first-order regression model to obtain data, and on this basis, then fit out the second-order regression model and obtain the final data. Or either one can be processed.

[0078] To ensure that one of the first-order / second-order regression models fitted is appropriate, it is necessary to use the regression analysis results to verify whether this regression model is significant, and the linear relationship strength between the independent variable and the dependent variable in the model can be known from the coefficient of determination. And the F-test of the lack of fit of the fitting is used to test whether the fitted regression model is appropriate.

[0079] Generally, when we are at a point far from the optimal condition on the response surface, there is only a small amount of curvature in the system and the first-order model is appropriate. At this time, we will quickly and effectively move towards the optimal point (optimum) along the improvement path. Once the best solution is reached and the existence of curvature is detected, a more refined second-order model is fitted.

[0080] The above first-order regression model and second-order regression model are linear regression algorithms.

[0081] In some preferred embodiments, as Figure 5 shown, step S400 includes:

[0082] S410. Construct a hybrid optimization strategy algorithm, and calculate a calculated value for obtaining the optimal solution of combined medication through this hybrid optimization strategy algorithm to obtain the maximum efficacy value;

[0083] Furthermore, the hybrid optimization strategy algorithm includes one of a genetic algorithm, a simulated annealing algorithm, and a particle swarm algorithm; the simulated annealing algorithm is used for global search, and the particle swarm algorithm is used for local optimization.

[0084] S420. Set the constraint conditions for the calculation method, and the constraint conditions include one or more of a single-dose range, a daily-dose range, and a safety limit.

[0085] In some preferred embodiments, as Figure 6 shown, step S100 includes:

[0086] S110. Construct a single-compartment PK calculation equation for a single drug.

[0087] Specifically in this step, a PK calculation equation needs to be constructed for each drug. For a single drug, the single-drug PK equation is constructed as follows:

[0088] Ci(t) = (Di * Ka) / (V * (Ka - Ke)) * (e^(-Ke * t) - e^(-Ka * t));

[0089] where Ci(t) is the blood drug concentration of drug i at time t;

[0090] Di is the administered dose;

[0091] Ka is the absorption rate constant;

[0092] Ke is the elimination rate constant;

[0093] V is the apparent volume of distribution;

[0094] t is the time after administration;

[0095] e is the base of the natural logarithm.

[0096] Furthermore, the two-compartment model PK equation is as follows:

[0097] Ci(t) = A 1e -αt + A 2e -βt ;

[0098] where:

[0099] A 1 =[D / (V 1 (α - β))]*[(k 21 - α)];

[0100] A 2 =[D / (V 1 (α - β))]*[(k 21 - β)];

[0101] α, β = 0.5{(k 12 + k 21 + k 10 ) ± [(k 12 + k 21 + k 10 )2 - 4k 21 k 10 )^(1 / 2)};

[0102] Parameter description:

[0103] C(t), drug concentration in the central compartment at time t;

[0104] D, administered dose;

[0105] V 1 , volume of distribution in the central compartment;

[0106] K 12 , first - order transfer rate constant from the central compartment to the peripheral compartment;

[0107] K 21 , first - order transfer rate constant from the peripheral compartment to the central compartment;

[0108] K 10 , first - order elimination rate constant for clearance from the central compartment;

[0109] α, rapid elimination rate constant in the distribution phase;

[0110] β, slow elimination rate constant in the elimination phase.

[0111] Furthermore, the PK equation for the three - compartment model:

[0112] C(t) = A 1e -αt + A 2e -βt + A 3e -γt ;

[0113] Where:

[0114] A 1 =[D / (V 1 )]*[(k 21-(α))(k 31 -(α)] / [(α - β)(α - γ)];

[0115] A 2 = [D / (V 1 )] * [(k 21 - β)(k 31 - β)] / [(β - α)(β - γ)];

[0116] A 3 = [D / (V 1 )] * [(k 21 - γ)(k 31 - γ)] / [(γ - α)(γ - β)].

[0117] α, β, γ are the three roots of the characteristic equation:

[0118] |λI - K| = 0;

[0119] where K is the rate constant matrix:

[0120] K = [-k 10 -k 12 -k 13 k 21 k 31

[0121] [k 12 -k 21 0]

[0122] [k 13 0 - k 31

[0123] Parameter description:

[0124] C(t): Drug concentration in the central compartment at time t;

[0125] D: Dose administered;

[0126] V 1 , Volume of distribution in the central compartment;

[0127] K 12 , Transfer rate constant from the central compartment to the first peripheral compartment;

[0128] K 21 , Transfer rate constant from the first peripheral compartment to the central compartment;

[0129] K 13 , Transfer rate constant from the central compartment to the second peripheral compartment;

[0130] K 31 , Transfer rate constant from the second peripheral compartment to the central compartment; ​​

[0131] K 10 , the elimination rate constant cleared from the central compartment;

[0132] α, the fastest elimination rate constant in the distribution phase;

[0133] β, the medium elimination rate constant in the distribution phase;

[0134] γ, the slowest elimination rate constant in the elimination phase;

[0135] S120. Construct a PK correction calculation equation for the single-drug PK calculation equation.

[0136] Specifically in this step, the PK calculation equations for each individual drug need to be corrected, and the correction calculation is as follows:

[0137] Ka_modified = Ka * (1 + ΣIij * Cj); Ke_modified = Ke * (1 + ΣIij * Cj);

[0138] Among them, Ka_modified is the corrected absorption rate constant considering drug interactions;

[0139] Ke_modified is the corrected elimination rate constant considering drug interactions;

[0140] Iij is the interaction coefficient of drug j on drug i;

[0141] Cj is the blood drug concentration of drug j.

[0142] S130. Calculate the AUC value based on the above single-drug PK calculation equation.

[0143] Specifically in this step, to obtain the AUC value, that is, to calculate the area under the drug concentration-time curve AUC, the calculation method is:

[0144] AUCi = ∫Ci(t)dt.

[0145] S140. Construct a response surface calculation model based on the above PK correction and AUC value.

[0146] Specifically in this step, according to the corrected blood drug concentration values and AUC calculation values in the above steps, a response surface model is constructed in the manner of the response surface method, and its calculation method is as follows:

[0147] y = β0 + Σβi * AUCi + Σβii * AUCi^2 + Σβij * AUCi * AUCj.

[0148] Among them, y is the pharmacodynamic response value;

[0149] β0 is the model intercept term;

[0150] βi is the first-order coefficient of AUCi;

[0151] βii: the second-order coefficient of AUCi;

[0152] βij is the interaction coefficient between AUCi and AUCj.

[0153] After the above steps, the PK equation, its numerical correction, and the calculation of AUC values are used as the factors for constructing the response surface method calculation, and the corresponding response surface model is constructed.

[0154] Furthermore, in some preferred embodiments, as Figure 7 shown, step S100 further includes designing an orthogonal test scheme, specifically:

[0155] S101. Determine the numerical values of drug dose levels.

[0156] In this step, the drug dose levels are delimited according to the available range obtained from the instructions or drug clinical trials. Assuming that the conventional treatment dose range of drug A is 50 mg to 200 mg / day, its dose levels can be divided into:

[0157] Low dose: 50 mg / day;

[0158] Medium dose: 100 mg to 150 mg / day;

[0159] High dose: 200 mg / day and above.

[0160] S102. Input the individualized parameter data of patients.

[0161] In this step, the individualized parameter data includes:

[0162] Basic physiological parameters: age (AGE), gender (SEX), weight (WT), body surface area (BSA), liver function indicators (ALT, AST, etc.) and kidney function indicators (SCr, CCr, etc.), etc.;

[0163] Genetic factors: CYP450 enzyme gene polymorphism, transporter gene polymorphism, and drug target genotype, etc.;

[0164] Disease status: degree of liver function impairment (Child-Pugh score, MELD score, ALBI score, etc.), degree of kidney function impairment (CKD-EPI formula, MDRD formula, KDIGO staging, RIFLE / AKIN criteria, etc.), and comorbidity situation (Charlson comorbidity index, Elixhauser comorbidity index, CIRS, etc.), etc.

[0165] S103. Establish an orthogonal test table.

[0166] Combining the above method steps, taking the combination of two drugs A and B in a one-compartment model as an example:

[0167] One-compartment model PK equation:

[0168] CA1(t) = (DAKaA) / (VA(KaA - KeA)) * (e^(-KeAt) - e^(-KaAt));

[0169] CB1(t) = (DBKaB) / (VB(KaB - KeB)) * (e^(-KeBt) - e^(-KaBt)).

[0170] Among them, CA(t) and CB(t) are the blood drug concentrations of drugs A and B at time t; DA and DB are the dosing doses of drugs A and B; KaA and KaB are the absorption rate constants of drugs A and B; KeA and KeB are the elimination rate constants of drugs A and B; VA and VB are the apparent distribution volumes of drugs A and B.

[0171] Correction considering interactions (taking the one-compartment model as an example):

[0172] KaA_modified(t) = KaA * (1 + IAB * CB1(t));

[0173] KeA_modified(t) = KeA * (1 + IAB * CB1(t));

[0174] KaB_modified(t) = KaB * (1 + IBA * CA1(t));

[0175] KeB_modified(t) = KeB * (1 + IBA * CA1(t)).

[0176] Among them, Ka_modified(t) and Ke_modified(t) represent the modified absorption and elimination rate constants at time t; since CB1(t) and CA1(t) are functions of time, the modified rate constants also change with time accordingly; this time-varying rate constant makes the pharmacokinetic equation become a nonlinear differential equation, and numerical methods need to be used to solve it; IAB is the interaction coefficient of drug B on drug A; IBA is the interaction coefficient of drug A on drug B.

[0177] The one-compartment model equations after correction:

[0178] dCA_modified / dt = KaA_modified(t)DA - KeA_modified(t)CA_modified;

[0179] dCB_modified / dt = KaB_modified(t)DB - KeB_modified(t)CB_modified.

[0180] Its solution is:

[0181] CA_modified(t) = (DAKaA_modified(t)) / (VA(KaA_modified(t) - KeA_modified(t)))*(e^(-KeA_modified(t)t) - e^(-KaA_modified(t)t));

[0182] CB_modified(t) = (DBKaB_modified(t)) / (VB(KaB_modified(t) - KeB_modified(t)))*(e^(-KeB_modified(t)t) - e^(-KaB_modified(t)t)).

[0183] Calculate the modified AUC:

[0184] AUCA = ∫CA_modified(t)dt;

[0185] AUCB = ∫CB_modified(t)dt.

[0186] Substitute into the response surface model:

[0187] y = β0 + βAAUCA + βBAUCB + βAAAUCA^2 + βBBAUCB^2 + βABAUCAAUCB.

[0188] AUCA and AUCB are functions of dose and individual pharmacokinetic parameters:

[0189] AUCA = f(DA, KaA, KeA, VA) = f(DA, θA);

[0190] AUCB = f(DB, KaB, KeB, VB) = f(DB, θB).

[0191] Where θA and θB are vectors of individual pharmacokinetic parameters, including: basic physiological parameters (age, gender, weight, body surface area, liver function indexes, kidney function indexes, etc.), genetic factors (CYP450 enzyme gene polymorphism, transporter gene polymorphism, drug target genotype, etc.), and pathological states (degree of liver insufficiency, degree of kidney insufficiency, and co - morbidity conditions). These individual factors affect PK parameters through the covariate model:

[0192] KaA = θKa,A * exp(ηKa,A + βage,KaAGE + βwt,KaWT + βgen,Ka * GEN);

[0193] KeA = θKe,A * exp(ηKe,A + βage,KeAGE + βwt,KeWT + βscr,Ke * SCR);

[0194] VA = θV,A * (WT / 70)^0.75 * exp(ηV,A).

[0195] Where: θ is the population parameter, η is the inter - individual variability, and β is the covariate effect coefficient

[0196] Therefore, the same dose may produce different AUCs in different individuals, and population differences need to be considered during optimization.

[0197] Steps for optimization solution:

[0198] Set the objective function to maximize the efficacy y

[0199] Set the constraints, which specifically include:

[0200] Single - dose range, Dmin ≤ Di,j ≤ Dmax (Dmin is the minimum allowable dose, Dmax is the maximum allowable dose, i represents the drug type, and j represents the j - th administration);

[0201] Daily - dose range, ΣDi,j ≤ DailyMax (DailyMax is the maximum daily dose);

[0202] Safety limit, Ci(t) ≤ Cmax (Cmax is the maximum safe concentration);

[0203] Dosing interval, Δti,j ≥ tmin,i (tmin,i is the minimum dosing interval time for drug i);

[0204] Course - of - treatment limit, ΣDi,j * N ≤ CourseMax (N is the number of days of the course of treatment, and CourseMax is the maximum course - of - treatment dose).

[0205] Furthermore, use numerical optimization methods to solve the optimal dose combination, including genetic algorithm solution, particle swarm optimization (PSO) solution, and simulated annealing algorithm (SA) solution.

[0206] Among them, the genetic algorithm solution process:

[0207] Encoding, encoding the doses DA and DB into binary chromosomes;

[0208] Initialization, randomly generate N individuals as the initial population;

[0209] Fitness evaluation, calculating the objective function value y of each individual;

[0210] Selection, using roulette wheel selection to choose excellent individuals;

[0211] Crossover, performing single-point crossover with a certain probability pc;

[0212] Mutation, performing gene mutation with a certain probability pm;

[0213] Repeat the iteration until the termination condition is met.

[0214] Among them, the particle swarm optimization (PSO) solution process:

[0215] Initialization, randomly generating the positions and velocities of N particles;

[0216] Calculating the fitness value of each particle;

[0217] Updating the individual best position pbest and the global best position gbest;

[0218] Updating the particle velocity according to the velocity update formula;

[0219] Updating the particle position according to the position update formula;

[0220] Repeat the iteration until convergence.

[0221] Among them, the simulated annealing algorithm (SA) solution process:

[0222] Initialization, randomly generating the initial solution and the initial temperature T

[0223] Generating a new solution at the current temperature;

[0224] Calculating the objective function difference ΔE;

[0225] Deciding whether to accept the new solution according to the Metropolis criterion;

[0226] Lowering the temperature T = αT (α is the temperature reduction coefficient);

[0227] Repeat the iteration until the temperature reaches the termination temperature.

[0228] Above, the genetic algorithm is suitable for large-scale search; the PSO algorithm converges quickly; the SA algorithm can jump out of the local optimum. Therefore, in specific implementation, first use SA for global search, and then use PSO for local optimization, which can take into account both the global search ability and the convergence speed. The genetic algorithm, the PSC algorithm, and the SA algorithm are all non-linear regression algorithms.

[0229] As described above, the advantage lies in using the AUC area of the drug metabolism process as the input variable, which can more accurately reflect the actual exposure of the drug in the body and the interaction effect; since AUC takes into account the whole process of drug absorption, distribution, metabolism and elimination, it can better characterize the actual bioavailability of the drug than simply using the administered dose, thus improving the accuracy and reliability of the combined drug dose optimization.

[0230] Combining the above methods and steps, the specific steps of the present invention for jointly incorporating the drug dose level and patient individual parameters into the orthogonal experimental design and introducing the AUC parameter into the response surface model are as follows:

[0231] Determination of factors and levels in the orthogonal experimental design, including: drug dose level (D) (drug A: DA1, DA2, DA3; drug B: DB1, DB2, DB3) and patient individual parameters (θ) (age AGE: young / middle-aged / elderly; weight WT: light / medium / heavy; CYP genotype: 1 / 1, 1 / 3, 3 / 3)

[0232] Establishing an individual PK model, including:

[0233] AUC expression of drug A:

[0234] AUCA = f(DA, θA) = f(DA, KaA, KeA, VA).

[0235] Among them, the relationship between individual parameters and covariates:

[0236] KaA = θka * (WT / 70)θwt * (AGE / 40)θage * θCYP;

[0237] KeA = θke * (WT / 70)θwt * (AGE / 40)θage * θCYP;

[0238] VA = θv * (WT / 70).

[0239] AUC expression of drug B:

[0240] AUCB = f(DB, θB) = f(DB, KaB, KeB, VB).

[0241] Orthogonal experimental design table (L18(37)):

[0242]

[0243] Construction of the response surface model, including:

[0244] Introducing AUC as a model parameter to obtain the calculation formula:

[0245] y = β0 + β1DA + β2DB + β3AGE + β4WT + β5CYP; β6AUCA + β7AUCB; β11DA2 + β22DB2; β12DADB + β67AUCAAUCB; β16DAAUCA + β27DBAUCB.

[0246] where AUCA = f(DA, θA); AUCB = f(DB, θB).

[0247] Model coefficient estimation includes:

[0248] Use the least squares method to fit a second-order regression model to the response quantity to obtain the linear relationship between the dose x and the response quantity;

[0249] Determine whether the second-order regression model is appropriate through the F-test;

[0250] If the model is appropriate, then solve the stationary point of the second-order regression model and determine whether the stationary point is a maximum point or a minimum point.

[0251] The optimization solution process is as follows:

[0252] Objective function: max(y);

[0253] Constraints: DA_min ≤ DA ≤ DA_max; DB_min ≤ DB ≤ DB_max; AUCA ≤ AUCA_max; AUCB ≤ AUCB_max.

[0254] The advantages of this solution: Consider the effects of both dose and individual parameters simultaneously; Reflect the actual exposure of the drug in the body through AUC; Establish a quantitative relationship between dose-exposure-effect; Achieve the individualized treatment goal of precise drug administration.

[0255] The combination of the above genetic algorithm and simulated annealing algorithm has more advantages. The specific process and analysis are as follows:

[0256] Use the genetic algorithm (GA) for preliminary search to determine the approximate range of the efficacy y.

[0257] Construct a fitness function based on the PK equation: f = w1*E(AUC) - w2*V(AUC) + w3*y;

[0258] Among them, E(AUC) is the expected drug exposure, representing the drug efficacy; V(AUC) is the inter-individual variability, representing the stability of the dosing regimen; y is a predefined clinical efficacy index (such as the symptom improvement rate); w1, w2, and w3 are weight coefficients to balance the drug efficacy, stability, and clinical efficacy; the chromosome encoding contains parameters such as [D1, D2, AGE, WT, GEN]; D1 and D2 are the doses of drugs A and B; AGE is the age group; WT is the weight group; GEN is the CYP genotype.

[0259] Advantages of the combination: Exploring a larger solution space through crossover mutation operations to find possible optimal efficacy intervals; Using genetic algorithms for global search can quickly locate possible optimal efficacy intervals, providing directions for subsequent optimization, and genetic algorithms can be computationally distributed, improving computational efficiency.

[0260] Taking experimental group 1 in the above orthogonal table as an example:

[0261] Initial chromosome: [DA1, DB1, elderly, heavy, *1 / *1]

[0262] Calculate the fitness (the following values are examples):

[0263] E(AUC) = Calculate the average AUC of the elderly-heavy weight-*1 / *1 genotype population = 42.3 mg·h / L;

[0264] V(AUC) = Calculate the variance of the AUC of this population (the population of the experimental group) = 6.8 (mg·h / L)^2;

[0265] y = Symptom improvement rate = 0.72; f = 0.4 * E(AUC) - 0.2 * V(AUC) + 0.4 * y = 15.64.

[0266] Crossover mutation: Cross with the chromosome [DA1, DB2, middle-aged, medium, *1 / *3] of experimental group 2 to generate a new solution

[0267] After 50 generations of evolution, assume that the global optimal solution is obtained: [DA1 = 100 mg, DB1 = 75 mg, elderly, heavy, *1 / *1].

[0268] Generally, the optimal solution under restricted conditions can be obtained. For example, assume the optimal solution after restricting the age group is as follows:

[0269] Optimal solution for the young group:

[0270] [DA1 = 85 mg, DB1 = 65 mg, young, medium, *1 / *1];

[0271] E(AUC) = 38.5 mg·h / L;

[0272] V(AUC) = 5.2 (mg·h / L)^2;

[0273] y = 0.78;

[0274] f = 0.4 * 38.5 - 0.2 * 5.2 + 0.4 * 0.78 = 14.76。

[0275] Optimal solution for the middle-aged group:

[0276] [DA1 = 95mg, DB1 = 70mg, middle-aged, medium, *1 / *1];

[0277] E(AUC) = 40.8 mg·h / L;

[0278] V(AUC) = 5.9 (mg·h / L)^2;

[0279] y = 0.75;

[0280] f = 0.4 * 40.8 - 0.2 * 5.9 + 0.4 * 0.75 = 15.34。

[0281] Optimal solution for the elderly group:

[0282] [DA1 = 75mg, DB1 = 55mg, elderly, medium, *1 / *1];

[0283] E(AUC) = 36.2 mg·h / L;

[0284] V(AUC) = 4.8 (mg·h / L)^2;

[0285] y = 0.71;

[0286] f = 0.4 * 36.2 - 0.2 * 4.8 + 0.4 * 0.71 = 13.85。

[0287] It can be seen from the above results that: the elderly group requires a lower dose to reduce the risk of adverse reactions; the middle-aged group can tolerate a higher dose to achieve the best therapeutic effect; the dose for the young group is in the middle, balancing the therapeutic effect and safety.

[0288] The simulated annealing algorithm (SA) is used for global optimization to optimize the dosing regimen while ensuring the therapeutic effect.

[0289] Based on multiple high-quality candidate solutions provided by GA as initial values, further optimization is carried out: GA provides multiple optimal solutions for different populations as the initial solution set for SA; each initial solution represents a local optimal region; SA can conduct a more detailed search around these regions.

[0290] Construct the state transition probability based on the GA solution:

[0291] P = exp(-ΔE / T), where ΔE = |AUCnew - AUCtarget| + λ|ynew - ytarget|;

[0292] Among them, P is the probability of accepting the new solution; T is the current temperature; AUCtarget is the target exposure; ytarget is the target efficacy value; λ is the balance coefficient.

[0293] Taking Experiment 3 in the above orthogonal table as an example:

[0294] Multiple candidate initial solutions provided by GA:

[0295] Solution 1, [DA2, DB2, young, light, *1 / *3], f = 16.80;

[0296] Solution 2, [DA2 + 0.2, DB2 - 0.1, young, light, *1 / *3], f = 16.75;

[0297] Solution 3, [DA2 - 0.1, DB2 + 0.15, young, light, *1 / *3], f = 16.72.

[0298] Select the optimal solution as the starting point of SA:

[0299] Initial solution: [DA2, DB2, young, light, *1 / *3]; E(AUC) = 45.6 mg·h / L; V(AUC) = 8.2 (mg·h / L)^2; y = 0.75; f = 0.4 * 45.6 - 0.2 * 8.2 + 0.4 * 0.75 = 16.80

[0300] The SA optimization process is complementary to the GA optimization. GA provides the results of large-scale search, and SA conducts fine search with small step sizes on this basis. The process is as follows:

[0301] Set the target values, AUCtarget = 50 mg·h / L; ytarget = 0.8;

[0302] Input the initial temperature, T0 = 100;

[0303] The first iteration, new solution, [DA2 + 0.1, DB2 - 0.05, young, light, *1 / *3], calculate ΔE = |47.2 - 50| + 0.5|0.77 - 0.8| = 3.95; P = exp(-3.95 / 100) = 0.961, accept the new solution;

[0304] Reduce the temperature, T = 100 * exp(-0.01 * 100 / 1000) = 99.0;

[0305] Continue the iteration.

[0306] This solution makes full use of the global search ability of GA and the local optimization characteristics of SA. The advantages are as follows: GA provides multiple potential high-quality solutions through population evolution; SA makes fine adjustments near these solutions to improve the optimization accuracy; the two algorithms complement each other's advantages, ensuring both the search scope and the optimization quality; this step further optimizes based on the GA solution through the SA algorithm, making the AUC and efficacy closer to the target values.

[0307] Introduce Bayesian Optimization (BO) for solution verification. As a verification step, predict the optimal efficacy. The process is as follows:

[0308] Construct a Gaussian process regression model, [y, AUC] ~ GP(μ(DA, DB, AGE, WT, CYP), k(DA, DB, AGE, WT, CYP)); μ(DA, DB, AGE, WT, CYP) is the mean function, describing the overall trend; k(DA, DB, AGE, WT, CYP) is the kernel function, describing the correlation structure;

[0309] Use the SA optimization result as prior knowledge, add the optimal solution obtained by SA and its neighborhood points to the training set, and update the hyperparameters (such as kernel function parameters) of the GP model;

[0310] The acquisition function guides the verification, EI(DA, DB, AGE, WT, CYP) = E[max(0, f(DA, DB, AGE, WT, CYP) - f(DA+, DB+, AGE, WT, CYP))]; f(DA+, DB+, AGE, WT, CYP) is the current optimal value; conduct key verification in the area with high EI values;

[0311] Iteratively verify and predict, sample new points in the candidate area; update the GP model; calculate the confidence interval of the optimal solution.

[0312] Regarding the problem of combination drug use, the existing technologies and their defects include: traditional empirical methods (lacking scientific quantitative basis), simple dose-effect models (ignoring the pharmacokinetic process and unable to reflect individual differences), statistical regression methods (insufficient individualization and unstable models when the sample size is insufficient), machine learning methods (strong data dependence and high overfitting risk in small sample cases), population pharmacokinetic methods (limited individual prediction accuracy and high requirements for sample representativeness), and Bayesian optimization methods (high computational complexity, relying on the quality of prior knowledge, and poor convergence in small samples).

[0313] Compared with the present invention, the defects existing in the prior art include:

[0314] 1. Only considering the static dose and ignoring the pharmacokinetic process (corresponding to the traditional empirical method and simple dose-effect model), it is thus impossible to accurately describe the processes of drug absorption, distribution, metabolism, and excretion in the body; lacking dynamic assessment of drug exposure; unable to reflect the impact of individual differences on drug metabolism; with low sample size requirements (20 - 30 cases), it is difficult to ensure the model accuracy.

[0315] 2. Unable to accurately evaluate drug interactions (corresponding to the statistical regression method), thus simply adding the effects of two drugs; only evaluating the drug efficacy based on the administered dose; unable to quantitatively describe the synergistic / antagonistic effects between drugs; a large number of samples (more than 100 cases) are required to establish a stable model; the surface fitting is inaccurate with small samples.

[0316] 3. Low degree of individualization and unified regimens (corresponding to machine learning methods, population pharmacokinetics methods, Bayesian optimization methods), thus not considering individual differences such as age, weight, genotype, etc.; lacking the mapping relationship between individual parameters and drug metabolism parameters; it is difficult to customize individualized dosing regimens; the sample size requirements are too high (200 - 500 cases), resulting in difficulties in clinical application; with a high risk of overfitting and poor model generalization ability in the case of small samples.

[0317] Combining the above embodiments, the beneficial effects of the present invention are as follows:

[0318] 1. The deep integration of the PK equation and the response surface model realizes the precise modeling of dynamic drug effects for the first time; compared with the prior art that only considers a single administered dose (such as: drug A, 100 mg twice a day, drug B, 50 mg twice a day), the present invention introduces the PK equation to describe the change in blood drug concentration when multiple drugs are used in combination (drug A: CA(t) = (DA * KaA) / (VA * (KaA - KeA)) * (e^(-KeA * t) - e^(-KaA * t)); drug B: CB(t) = (DB * KaB) / (VB * (KaB - KeB)) * (e^(-KeB * t) - e^(-KaB * t)); interaction: KeA = KeA0 * (1 + α * CB(t)), where α is the inhibition coefficient of drug B on the metabolism of A). Specifically, for example, when the antihypertensive drug A (CYP3A4 substrate) is used in combination with the antifungal drug B (a potent CYP3A4 inhibitor), the traditional fixed-dose regimen cannot evaluate the increase in blood drug concentration caused by drug interactions, while the present invention can predict that the AUC of drug A increases by 2 - 3 times through the PK equation. Accordingly, the dose of drug A is reduced to 30 - 50% of the original dose (twice a day, 30 - 50 mg each time), and at the same time, the dosing regimen of drug B is adjusted to avoid adverse reactions.

[0319] 2. Significantly improved the scientific nature of combined drug use and achieved precise optimization. Compared with the prior art which simply superimposes the effects of two drugs and only evaluates the drug efficacy based on the administered dose; only models based on drug doses DA and DB cannot reflect the actual exposure of the drug in the body; it is more reasonable to evaluate the drug efficacy through the AUC index in the present invention because AUC reflects the overall exposure of the drug in the body, while the dose only represents the administered amount; the same dose may result in different AUCs among different individuals, leading to differences in drug efficacy; AUC takes into account the whole process of drug absorption, distribution, metabolism, and clearance; the drug interaction can be quantitatively described through the "AUCA*AUCB" interaction term; by introducing the AUC index, first-order and second-order response surface models are established. First-order model: y = β0 + β1AUCA + β2AUCB + β12AUCA*AUCB; Second-order model: y = β0 + β1AUCA + β2AUCB + β12AUCA*AUCB + β11AUCA2 + β22AUCB2, where AUCA and AUCB are calculated through the PK equation. Furthermore, in specific implementations, for example, in the case of co-administering a CYP enzyme inhibitor, the traditional method based on the dose cannot predict drug interactions, while the present invention can accurately reflect the increase in drug exposure caused by enzyme inhibition through AUC, thereby timely adjusting the dosing regimen; for the combination of antihypertensive drugs A and B, the traditional method only considers the administered dose and cannot accurately reflect the change in blood drug concentration caused by individual differences, while the present invention can precisely quantify the actual exposure of the drug in the body through the AUC index, thereby finding the optimal individualized dosing regimen.

[0320] 3. The innovative modeling method based on AUC overcomes the limitations of simple dose indicators. Compared with the prior art which only considers the administered dose and cannot reflect the actual exposure of the drug in the body, the present invention introduces the AUC index to evaluate drug exposure, and the correlation between AUC and efficacy is stronger; AUC directly reflects the exposure of the drug at the target site and is closely related to the pharmacodynamic effect; the same dose may result in different AUCs, leading to differences in efficacy; AUC takes into account the influence of individual differences on drug absorption, distribution, metabolism, and clearance. In specific implementations, for example, for patients with liver insufficiency, the traditional fixed-dose regimen may lead to drug accumulation, while in the present invention, when it is found through monitoring that AUC increases by more than 50%, the dose is timely adjusted to 50 - 75% of the normal dose to avoid adverse reactions.

[0321] 4. Provide a more comprehensive assessment of the drug exposure-effect relationship. Compared with the prior art that only focuses on peak concentration and trough concentration, which has a poor correlation with efficacy and is difficult to evaluate the relationship between overall drug exposure and efficacy, the present invention establishes a more scientific exposure-effect relationship through AUC. AUC has a better correlation with clinical efficacy, can accurately reflect the therapeutic equivalence of different dosing regimens, and individualized dosing based on AUC is more likely to achieve the target efficacy. Establish a quantitative relationship between AUC and clinical endpoints to determine the AUC range corresponding to the optimal therapeutic effect. Establish a dose-effect curve of the adverse reaction risk and AUC. In specific implementation, for example, in antibiotic treatment, the traditional method only monitors the ratio of peak and trough concentrations to MIC, with a poor correlation, while the present invention can more accurately predict clinical efficacy and optimize the dosing regimen through the AUC / MIC ratio.

[0322] 5. Innovative individualized orthogonal experimental design, and systematically considering multi-dimensional individual characteristics. Compared with the prior art that only considers 1-2 main factors, such as age and weight, the present invention uses the L18(3^5) orthogonal array to consider 5 key factors simultaneously (L18 represents a Latin square orthogonal array containing 18 experiments, and each of the 5 factors has 3 levels). In specific implementation, for example, in anticoagulant drugs, the traditional method only adjusts the dose according to weight, while the present invention simultaneously considers factors such as age, liver and kidney function, and gene polymorphism, significantly reducing the bleeding risk.

[0323] 6. Optimize the experimental efficiency and reduce the sample size collected. Compared with the completely randomized experimental design in the prior art, which requires a large number of samples, resulting in high experimental costs, difficult actual operation, and low feasibility, the present invention optimizes the strategy through orthogonal experiments, has a uniformly distributed level combination, minimizes the number of experiments, and has the advantage of ensuring the representativeness of the experiment. And through the stratified sampling scheme, stratify according to key characteristics to ensure the sample size of each subgroup. In specific implementation, compared with the random experiments in the prior art that require hundreds of samples, the present invention only needs 18 groups of experiments through the L18(3^5) orthogonal array to obtain reliable conclusions.

[0324] As Figure 8 shown, a device 1000 is used to execute the method for determining the combination medication described in any of the above embodiments, and includes: a model architecture module 1001 for constructing a computational model, an experimental design module 1002 for constructing an experimental scheme, an optimization solving module 1003 for calculating the optimal solution of the combination medication, and a distributed computing module 1004 for accelerating the calculation process.

[0325] A readable storage medium stores a program or instruction thereon, and when the program or instruction is executed by a processor, the steps of the method in any one of the embodiments are implemented.

[0326] The preferred specific embodiments of the present invention have been described in detail above. It should be understood that those of ordinary skill in the art can make many modifications and variations based on the concept of the present invention without creative efforts. Therefore, all technical solutions that can be obtained by those skilled in the art in the technical field according to the concept of the present invention through logical analysis, reasoning or limited experiments on the basis of the prior art should fall within the protection scope determined by the claims.

Claims

1. A method for determining a drug combination, characterized in that: include: Constructing a computational model, wherein the computational model is constructed according to response surface methodology; Obtaining a reaction value, and inputting the reaction value into the calculation model; A linear regression algorithm or / and a nonlinear regression algorithm suitable for the relationship between the response variable and the response value; The optimal solution for drug combination is obtained, and the optimal solution is obtained by using the linear regression algorithm and / or the nonlinear regression algorithm.

2. A method for determining a drug combination according to claim 1, characterized in that: The steps are adapted to a linear regression algorithm or / and a nonlinear regression algorithm for the relationship between the response variable and the response value, and include: Fitting a first-order regression model and / or a second-order regression model; Inputting the response value into the first-order regression model or / and the second-order regression model, and obtaining the variable value; A lack of suitability detection algorithm adapted to the first-order regression model and / or the second-order regression model; Inputting the variable value into the deficiency test algorithm and obtaining the test result; Inputting the verification result, and judging the adaptability of the first-order regression model and / or the second-order regression model according to the verification result; Obtain a first-order regression model or / and a second-order regression model that meets adaptability.

3. A method for determining a drug combination as claimed in claim 2, characterized in that: include: Calculating a stable point value, wherein the stable point value is obtained by calculating the first-order regression model or / and the second-order regression model that satisfies adaptability; It is determined whether the stable point value satisfies the optimal solution for the combined medication. If so, the corresponding stable point value is adopted.

4. The method for determining a drug combination according to claim 1, wherein: The steps for obtaining the best solution for drug combination therapy include: Constructing a hybrid optimization strategy algorithm, and calculating a calculated value for obtaining an optimal solution for combined drug therapy through the hybrid optimization strategy algorithm to obtain a maximum therapeutic effect value; Constraints for the calculation method are set, and the constraints include one or more of a single dose range, a daily dose range, and a safety limit.

5. A method for determining a drug combination according to claim 4, characterized in that: The hybrid optimization strategy algorithm includes one of a genetic algorithm, a simulated annealing algorithm and a particle swarm algorithm, and the simulated annealing algorithm is used for global search, and the particle swarm algorithm is used for local optimization.

6. The method for determining a drug combination according to claim 1, wherein: The step constructs a computational model, and the computational model is constructed according to the response surface methodology and includes: Construct single-drug PK calculation equation; Constructing a PK correction calculation equation for the single-drug PK calculation equation; Calculate the AUC value based on the above single-drug PK calculation equation; A response surface calculation model based on the above PK corrections and AUC values ​​was constructed.

7. A method for determining a drug combination according to claim 6, characterized in that: The step constructs a computational model, and the computational model is constructed according to the response surface methodology and includes: Determine numerical values ​​for drug dosage levels; Enter patient-specific parameter data; Create an orthogonal test table.

8. A method for determining a drug combination according to claim 7, characterized in that: The individualized parameter data include: basic physiological parameters, genetic factors and disease status.

9. A device, characterized in that: The method for determining the drug combination described in any one of claims 1 to 8 above includes: a model construction module for constructing a calculation model, an experimental design module for constructing an experimental plan, an optimization solution module for calculating the best solution for the drug combination, and a distributed computing module for accelerating the calculation process.

10. A readable storage medium storing a program or instruction, wherein the program or instruction, when executed by a processor, implements the steps of the method of any one of claims 1 to 8.

Images