如何使用耦合方程拟合非线性回归

问：

我想用一些实验数据来拟合这个方程：nls

Model <- as.formula (F ~  (dhg*K1*G+dhg2*K1*K2*G^2)/(1 + K1*G+K1*K2*G^2)  - d0)

但已从以下位置获得：G

K1*K2*G^3*+K1*(2*K2*Ht-K2*Gt.M+1)*G^2+K1*(Ht-Gt.M)+1*G-Gt = 0

有没有办法将函数包含在公式中（因此可以在每次迭代中求解，例如，通过使用 polyroot 获取根）？例如G

G <- function(K1, K2, Ht, Gt){
    D <- K1*K2
    C <- K1*(2*K2*Ht-K2*Gt+1)
    B <- K1*(Ht-Gt)+1
    A <- -Gt
    G <- Re(polyroot(z=c(A,B,C,D)))[which(Re(polyroot(z=c(A,B,C,D))) > 0)]
    return(G)
}

我没有设法在公式中包含 G

编辑

感谢 Onyambu 的回答，我制作了以下工作示例：

Data <- data.frame(
    Gt.M =c(0.000000000,0.001085104,0.002170209,0.003255313,0.004340418,0.005425522,0.006510627,0.007595731,0.008680836,0.009765940,0.010851045,0.011936149,0.013021254,0.014106358,0.015191463,0.016276567),
    Fitted_d.max =c(7.167307, 7.162498, 7.158839, 7.155603, 7.152150, 7.149008, 7.145748, 7.142980,7.139988, 7.137311, 7.134371, 7.131921, 7.129423, 7.127004, 7.124224, 7.121391))

G <- function(K1, K2, Ht, Gt.M){
    D <- K1*K2
    C <- K1*(2*K2*Ht-K2*Gt.M+1)
    B <- K1*(Ht-Gt.M)+1
    A <- -Gt.M
    G <- Re(polyroot(z=c(A,B,C,D)))[which(Re(polyroot(z=c(A,B,C,D))) > 0)]
    return(G)
}

Model2_1 <- as.formula (Fitted_d.max ~ d0 + (dhg*K1*G(K1, K2, Gt.M,Ht=1E-3)+dhg2*K1*K2*G(K1, K2, Gt.M,Ht=1E-3)^2)/(1 + K1*G(K1, K2, Gt.M,Ht=1E-3)+K1*K2*G(K1, K2, Gt.M,Ht=1E-3)^2))

nls(Model2_1,
        data = Data, 
        start = list(K1 = 1, K2=1, dhg=0.1, dhg2=0.1,d0=7.3))

但是，它会返回错误

Error in nls(Model2_1, data = Data, start = list(K1 = 1, K2 = 1, dhg = 0.1,  : 
  step factor 0.000488281 reduced below 'minFactor' of 0.000976562
In addition: There were 50 or more warnings (use warnings() to see the first 50)

并使用 nlstools 包中的 preview（） 函数，它返回：

Error in data.frame(..., check.names = FALSE) : 
  arguments imply differing number of rows: 16, 23

恐怕我不确定问题出在哪里

谢谢！

编辑2

我通过改变计算方式来让它工作G

library('VGAMextra')

Ft <- function(K1, K2, Gt.M,Ht,d0,dhg,dhg2){
G <- function(K1, K2, Gt.M,Ht){
F <- function(K1, K2, Gt.M,G,Ht) K1*K2*G^3 + K1*(2*K2*Ht-K2*Gt.M+1)*G^2 + (K1*(Ht-Gt.M)+1)* G -Gt.M
Fdev <- function(K1, K2, Gt.M,G,Ht) 3*K1*K2*G^2 + 2*K1*(2*K2*Ht-K2*Gt.M+1)*G + (K1*(Ht-Gt.M)+1)
newtonRaphson.basic(f = F, fprime  = Fdev, a = 0, b= Gt.M, K1=K1, K2=K2, Gt.M=Gt.M,Ht=Ht,nmax=100)}
Gc <- G(K1=K1, K2=K2, Gt.M = Gt.M, Ht=Ht)
Fitted_d.max <- d0 + (dhg*K1*Gc+dhg2*K1*K2*Gc^2)/(1 + K1*Gc+K1*K2*Gc^2)
}

 nls(Fitted_d.max~Ft(K1, K2, Ht=0.5E-3, Gt.M, d0 = Data[1,2],dhg,dhg2),
        data = Data, 
        start = list(K1 = 1500, K2=5, dhg=-0.1, dhg2=-0.19),
        lower= c(0,0, -1,-1),
        trace=TRUE, algorithm="port")

似乎在做这件事。