Negative exponentielle Anpassung: Kurve sieht zu hoch aus

Question

Ich versuche, ein negatives Exponential zu einigen Daten in R zu passen, aber die angepasste Linie sieht zu hoch im Vergleich zu den Daten, während die Passform mit Excel eingebauten Power Fit aussieht glaubwürdiger. Kann mir jemand sagen warum? Ich habe versucht, mit der nls() Funktion und auch optim() und ähnliche Parameter von beiden dieser Methoden, aber die passt für beide hoch.

x <- c(5.96, 12.86, 8.40, 2.03, 12.84, 21.44, 21.45, 19.97, 8.92, 25.00, 19.90, 20.00, 20.70, 16.68, 14.90, 26.00, 22.00, 22.00, 10.00, 5.70, 5.40, 3.20, 7.60, 0.59, 0.14, 0.85, 9.20, 0.79, 1.40, 2.68, 1.91) 
    y <- c(5.35, 2.38, 1.77, 1.87, 1.47, 3.27, 2.01, 0.52, 2.72, 0.85, 1.60, 1.37, 1.48, 0.39, 2.39, 1.83, 0.71, 1.24, 3.14, 2.16, 2.22, 11.50, 8.32, 38.98, 16.78, 32.66, 3.89, 1.89, 8.71, 9.74, 23.14) 

    xy.frame <- data.frame(x,y) 

    nl.fit <- nls(formula=(y ~ a * x^b), data=xy.frame, start = c(a=10, b=-0.7)) 

    a.est <- coef(nl.fit)[1] 
    b.est <- coef(nl.fit)[2] 

    plot(x=xy.frame$x,y=xy.frame$y) 

    # curve looks too high 
    curve(a.est * x^b.est , add=T) 
    # these parameters from Excel seem to fit better 
    curve(10.495 * x^-0.655, add=T) 

# alternatively use optim() 
    theta.init <- c(1000,-0.5, 50) 

    exp.nll <- function(theta, data){ 
     a <- theta[1] 
     b <- theta[2] 
     sigma <- theta[3] 
     obs.y <- data$y 
     x <- data$x 
     pred.y <- a*x^b 
     nll <- -sum(dnorm(x=obs.y, mean=pred.y , sd=sigma, log=T)) 
     nll 
    } 

    fit.optim <- optim(par=theta.init,fn=exp.nll,method="BFGS",data=xy.frame) 

    plot(x=xy.frame$x,y=xy.frame$y) 

    # still looks too high 
    curve(a.est * x^b.est, add=T) 

Answer 1

Der Grund, warum Sie das unerwartete Verhalten zu sehen sind, ist, dass die Kurven, die „zu hoch“ tatsächlich aussehen viel geringere Summen der quadrierten Fehler als die Kurven von Excel:

# Fit from nls 
sum((y - a.est*x^b.est)^2) 
# [1] 1588.313 

# Fit from excel 
sum((y - 10.495*x^ -0.655)^2) 
# [1] 1981.561 

Der Grund nls fa Vor der höheren Kurve ist es, dass große Fehler bei kleinen x-Werten auf Kosten von etwas größeren Fehlern mit großen x-Werten vermieden werden. Eine Möglichkeit, dies zu adressieren könnte sein, eine Protokoll-Log-Transformation anzuwenden:

mod <- lm(log(y)~log(x)) 
(a.est2 <- exp(coef(mod)["(Intercept)"])) 
# (Intercept) 
# 10.45614 
(b.est2 <- coef(mod)["log(x)"]) 
#  log(x) 
# -0.6529741 

Diese sind sehr nahe an die Koeffizienten von Excel und ergibt eine optisch ansprechende Passform (trotz der schlechteren Leistung auf dem Summe-von- eckige Fehler metrisch):