I'm trying to write functions that wrap basic gamlss modeling (because I'm preparing a higher-level, application-specific package). But I've run into issues wrapping predict.
The example from the predict.gamlss
documentation works well:
suppressMessages(library(gamlss))
data(aids)
a<-gamlss(y~poly(x,3)+qrt, family=PO, data=aids)
#> GAMLSS-RS iteration 1: Global Deviance = 416.8014
#> GAMLSS-RS iteration 2: Global Deviance = 416.8014
newaids<-data.frame(x=c(45,46,47), qrt=c(2,3,4))
predict(a, newdata = newaids)
#> Warning in predict.gamlss(a, newdata = newaids): There is a discrepancy between the original and the re-fit
#> used to achieve 'safe' predictions
#>
#> [1] 6.021313 6.232357 6.155965
But as soon as I wrap it, the code fails:
gamlss_fit <- function(input) {
gamlss(y~poly(x,3)+qrt, family=PO, data=input)
}
gamlss_predict <- function(fit, input) {
predict(fit, newdata = input)
}
fit <- gamlss_fit(aids)
#> GAMLSS-RS iteration 1: Global Deviance = 416.8014
#> GAMLSS-RS iteration 2: Global Deviance = 416.8014
gamlss_predict(fit, newaids)
#> Error in eval(Call$data): object 'input' not found
I vaguely suspect there's an issue managing the assumed environment
. I've tried looking at the implementation of predict.gamlss
, but I haven't been able to follow it. The problem seems related to the question here, but the original newdata
call works in my case. I'd appreciate any suggestions.
I was able to solve the same issue by add the input dataset to the global environment within the function. Try:
gamlss_predict <- function(fit, input) {
input <<- input
predict(fit, newdata = input)
rm(input)
}