Vorhersage mit zufälligen Effekten in mgcv gam

Ich bin daran interessiert, den gesamten Fischfang mit gam in mgcv zu modellieren, um einfache zufällige Effekte für einzelne Schiffe zu modellieren (die im Laufe der Zeit wiederholte Fahrten in der Fischerei unternehmen). Ich habe 98 Probanden, also dachte ich, ich würde Gam anstelle von Gamm verwenden, um die zufälligen Effekte zu modellieren. Mein Modell ist:

modelGOM <- gam(TotalFish ~ factor(SetYear) + factor(SetMonth) + factor(TimePeriod) +     
s(SST) + s(VesselID, bs = "re", by = dum) + s(Distance, by = TimePeriod) + 
offset(log(HooksSet)), data = GOM, family = tw(), method = "REML")

Ich habe den Zufallseffekt mit bs = "re" und by = dum codiert (ich habe gelesen, dass ich damit die Gefäßeffekte bei ihren vorhergesagten Werten oder Null vorhersagen kann). "dum" ist ein Vektor von 1.

Das Modell läuft, aber ich habe Probleme bei der Vorhersage. Ich habe eines der Schiffe für die Vorhersagen (Schiff 21) und Durchschnittswerte für alles andere außer dem für Vorhersagen interessanten Prädiktor (Entfernung) ausgewählt.

data.frame("Distance"=seq(min(GOM$Distance),max(GOM$Distance),length = 100),
                             "SetYear" = '2006',
                             "SetMonth" = '6',
                             "TimePeriod" = 'A',
                             "SST" = mean(GOM$SST),
                             "VesselID" = 'Vessel21', 
                             "dum" = '0', #to predict without vessel effect
                             "HooksSet" = mean(GOM$HooksSet))

pred_GOM_A_Swordfish <- predict(modelGOM, grid.bin.GOM_A_Swordfish, type = "response", 
se = T)

Der Fehler, den ich bekomme, ist:

Error in Predict.matrix.tprs.smooth(object, dk$data) : 
    NA/NaN/Inf in foreign function call (arg 1)
    In addition: Warning message:
    In Ops.factor(xx, object$shift[i]) : - not meaningful for factors

Ich denke, das wird aufgerufen, weil VesselID ein Faktor ist, aber ich benutze es glatt für die zufälligen Effekte.

Ich konnte die Verwendung von Gam ohne die einfachen zufälligen Effekte (bs = "re") erfolgreich vorhersagen.

Können Sie Ratschläge geben, wie Sie dieses Modell ohne den VesselID-Begriff vorhersagen können (und es dennoch in die Anpassung einbeziehen)?

Vielen Dank!

Ab Version 1.8.8 von mgcv predict.gam wurde ein excludeArgument gewonnen, das das Nullstellen von Begriffen im Modell, einschließlich zufälliger Effekte, bei der Vorhersage ohne den zuvor vorgeschlagenen Dummy-Trick ermöglicht.

• predict.gamund predict.bamakzeptieren Sie nun ein 'exclude'Argument, mit dem Begriffe (z. B. zufällige Effekte) zur Vorhersage auf Null gesetzt werden können. Aus Effizienzgründen werden glatte Begriffe, die nicht in termsoder in excludesind, nicht mehr ausgewertet und stattdessen auf Null gesetzt oder nicht zurückgegeben. Siehe ?predict.gam.

library("mgcv")
require("nlme")
dum <- rep(1,18)
b1 <- gam(travel ~ s(Rail, bs="re", by=dum), data=Rail, method="REML")
b2 <- gam(travel ~ s(Rail, bs="re"), data=Rail, method="REML")

head(predict(b1, newdata = cbind(Rail, dum = dum)))    # ranefs on
head(predict(b1, newdata = cbind(Rail, dum = 0)))      # ranefs off
head(predict(b2, newdata = Rail, exclude = "s(Rail)")) # ranefs off, no dummy

> head(predict(b1, newdata = cbind(Rail, dum = dum)))    # ranefs on
       1        2        3        4        5        6 
54.10852 54.10852 54.10852 31.96909 31.96909 31.96909  
> head(predict(b1, newdata = cbind(Rail, dum = 0)))      # ranefs off
   1    2    3    4    5    6 
66.5 66.5 66.5 66.5 66.5 66.5
> head(predict(b2, newdata = Rail, exclude = "s(Rail)")) # ranefs off, no dummy
   1    2    3    4    5    6 
66.5 66.5 66.5 66.5 66.5 66.5

Älterer Ansatz

Simon Wood hat das folgende einfache Beispiel verwendet, um zu überprüfen, ob dies funktioniert:

library("mgcv")
require("nlme")
dum <- rep(1,18)
b <- gam(travel ~ s(Rail, bs="re", by=dum), data=Rail, method="REML")
predict(b, newdata=data.frame(Rail="1", dum=0)) ## r.e. "turned off"
predict(b, newdata=data.frame(Rail="1", dum=1)) ## prediction with r.e

Welches funktioniert bei mir. Gleichfalls:

dum <- rep(1, NROW(na.omit(Orthodont)))
m <- gam(distance ~ s(age, bs = "re", by = dum) + Sex, data = Orthodont)
predict(m, data.frame(age = 8, Sex = "Female", dum = 1))
predict(m, data.frame(age = 8, Sex = "Female", dum = 0))

funktioniert auch.

Ich würde also überprüfen, ob die von Ihnen bereitgestellten Daten Ihrer newdataMeinung nach das Problem sind, mit dem Sie möglicherweise nicht zu VesselIDtun haben. Der Fehler stammt von der Funktion, die von den predict()Aufrufen in den obigen Beispielen aufgerufen worden wäre , und Rail ist ein Faktor in der erstes Beispiel.