R Anpassungsgleichung (Nikol'sky-Eisenman-Gleichung), um Daten

Question

verwendet I SigmaPlot verwenden, um eine modifizierte Version der Nikol'sky-Eisenman-Gleichung der Form Fitfunktion

y = P1 + P2 * log(10^(-x) + P3) 

unter Verwendung der globalen Kurve zu passen. Details zu den Parametern finden Sie im Sigmaplot-Bericht unten. Ich möchte dies in R. jetzt tun

Einige Daten:

pNO3 <- c(1.1203, 2.0410, 3.0155, 4.0048, 4.3045, 5.0, 6.0) 
mV <- c(45.2, 100.9, 160.9, 215.7, 231.5, 244.5, 257.4) 
data <- data.frame(pNO3, mV) 
plot(data$pNO3, data$mV) 

Die Grafik und Bericht von SigmaPlot für die obigen Daten erzeugt wird unten gezeigt. Kann mir jemand in die richtige Richtung zeigen, wie man etwas ähnliches in R erzeugt?

NonLinear Regression - Global Curve Fitting  Wednesday, May 01, 2013, 13:04:55 
Data Source: Data 1 in Notebook1 
Equation: User-Defined, Nicolsky Eisenman 
f=P1+P2*log(10^(-x)+P3) 
Data Set Specifications: 
Data Set Independent Variable Dependent Variable 
1   Column 3  Column 7 
Global Parameters: 
A Global Parameter is shared across all data sets. 
Global Goodness of Fit: 
R Rsqr Adj Rsqr Standard Error of Estimate 
0.9997 0.9994 0.9991  2.4421 
Analysis of Variance: 
Analysis of Variance: 
    DF SS MS 
Regression 3 264242.5551 88080.8517 
Residual 4 23.8549 5.9637 
Total 7 264266.4100 37752.3443 
Corrected for the mean of the observations: 
    DF SS MS F P 
Regression 2 38844.3822 19422.1911 3256.7192 <0.0001 
Residual 4 23.8549 5.9637 
Total 6 38868.2371 6478.0395 
Statistical Tests: 
Normality Test (Shapiro-Wilk)   Passed (P = 0.4003) 
W Statistic= 0.9106 Significance Level = 0.0500 
Constant Variance Test  Passed (P = 0.1209) 
Number of Observations = 7 
Rsqr = 0.9994 
Residual Sum of Squares = 23.8549 
Parameter Estimates: 
    Coefficient Std. Error t P 
P1 -24.3265  3.3330 -7.2987 0.0019 
P2 -61.7088  1.2861 -47.9796 <0.0001 
P3 2.8351E-005 4.6040E-006 6.1579 0.0035 
Fit Equation Description: 
[Variables] 
f0_x = col(3) 
f0_y = col(7) 
[Parameters] 
f0_P1 = 0 ' {{previous: -24.3265}} 
f0_P2 = -5 ' {{previous: -61.7088}} 
f0_P3 = 0 ' {{previous: 2.8351e-005}} 
[Equation] 
f0 = f0_P1+f0_P2*log(10^(-f0_x)+f0_P3) 
fit f0 to f0_y 
[Constraints] 
[Options] 
tolerance=0.000100 
stepsize=100 
iterations=100 
Number of Iterations Performed = 4 

Answer 1

Ihr Passkriterium Unter der Annahme, das summierte quadratische Fehler zu minimieren, Sie nls verwenden können, aber Sie tun, einen fairen Startwert benötigen. Ich habe keine Ahnung, was für Ihre Parameter vernünftig ist, also habe ich eine Weile gebraucht, bis ich die Parameter aus Ihrem sigmaplot-Beispiel kopiert habe, was meiner Meinung nach für den Datensatz angemessen war, der diesem Datensatz ähnlich sein könnte. Wie auch immer, wenn Sie wissen, was die Parameter bedeuten, können Sie wahrscheinlich vernünftige Startwerte schätzen.

> start=list(P1=-24,P2=-61,P3=2.8e-5) 
> m = nls(formula= mV ~ P1 + P2 * log(10^(-pNO3) + P3),data=data,start=start) 
> summary(m) 

Formula: mV ~ P1 + P2 * log(10^(-pNO3) + P3) 

Parameters: 
    Estimate Std. Error t value Pr(>|t|)  
P1 -1.420e+01 4.642e+00 -3.059 0.055 . 
P2 -2.732e+01 9.257e-01 -29.514 8.54e-05 *** 
P3 8.417e-05 1.818e-05 4.630 0.019 * 

Sie die Daten dargestellt werden können, und die glatte Kurvenanpassung durch einen neuen Satz von pNO3 Maßnahmen zu schaffen:

plot(data$pNO3,data$mV) 
newdata = data.frame(pNO3=seq(1,6,len=100)) 
lines(newdata$pNO3,predict(m, newdata=newdata)) 

Beachten Sie, dass "log" in R natürlichen Logarithmus ist, wenn Sie wollen Log auf Basis 10 dann verwenden Sie log10 - dies ändert P2 ein bisschen, um etwa -62 statt -27 wie oben ...

Mit y unsere neuen Daten und mit „log10“ anstelle von „log“ in der Formel Ausdruck I erhalten:

> m10 = nls(formula= mV ~ P1 + P2 * log10(10^(-pNO3) + P3),data=data,start=start) 
> summary(m10) 

Formula: mV ~ P1 + P2 * log10(10^(-pNO3) + P3) 

Parameters: 
    Estimate Std. Error t value Pr(>|t|)  
P1 -2.433e+01 3.334e+00 -7.298 0.00187 ** 
P2 -6.171e+01 1.286e+00 -47.972 1.13e-06 *** 
P3 2.835e-05 4.605e-06 6.157 0.00353 ** 

, die wie Ihre Sigmaplot Ausgabe sieht:

Parameter Estimates: 
    Coefficient Std. Error t P 
P1 -24.3265  3.3330 -7.2987 0.0019 
P2 -61.7088  1.2861 -47.9796 <0.0001 
P3 2.8351E-005 4.6040E-006 6.1579 0.0035 

Answer 2

Mit dem "plinear" Algorithmus von nls nur die Parameter die nichtlinear eingegeben werden müssen Startwerte. Beachten Sie, dass plinear eine Änderung in der Formel erfordert, wie gezeigt, dieses Modell zu repräsentieren:

fo <- mV ~ cbind(log10(10^(-pNO3) + P3), 1) 
fm <- nls(fo, data, start = c(P3 = 0), alg = "plinear") 
summary(fm) 

geben:

Formula: mV ~ cbind(log10(10^(-pNO3) + P3), 1) 

Parameters: 
     Estimate Std. Error t value Pr(>|t|)  
P3  2.84e-05 4.60e-06 6.16 0.0035 ** 
.lin1 -6.17e+01 1.29e+00 -47.97 1.1e-06 *** 
.lin2 -2.43e+01 3.33e+00 -7.30 0.0019 ** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Residual standard error: 2.44 on 4 degrees of freedom 

Number of iterations to convergence: 7 
Achieved convergence tolerance: 2.94e-06 

Und wir können es so Grundstück:

plot(data) 
lines(fitted(fm) ~ pNO3, data, col = "red")