GeNetIt (version 0.1-6)

predict.gravity: Predict gravity model

Description

predict method for class "gravity"

Usage

# S3 method for gravity
predict(
  object,
  newdata,
  groups = NULL,
  back.transform = c("none", "simple", "Miller", "Naihua"),
  ...
)

Value

Vector of model predictions

Arguments

object

Object of class gravity

newdata

New data used for obtaining the predictions, can be a data.frame or nffGroupedData

groups

Grouping factor acting as random effect. If used, must match levels used in model, otherwise leave it null and do not convert to groupedData

back.transform

Method to back transform data, default is none and log predictions will be returned.

...

Arguments passed to predict.lme or predict.lm

Author

Jeffrey S. Evans <[e-mail]> and Melanie A. Murphy <[e-mail]>

Details

Please note that the entire gravity equation is log transformed so, your parameter space is on a log scale, not just y. This means that for a meaningful prediction the "newdata" also needs to be on a log scale.

For the back.transform argument, the simple back-transform method uses the form exp(y-hat)0.5*variance whereas Miller uses exp(sigma)*0.5 as the multiplicative bias factor. Naihua regresses y~exp(y-hat) with no intercept and uses the resulting coefficient as the multiplicative bias factor. The Naihua method is intended for results with non-normal errors. You can check the functional form by simply plotting y (non-transformed) against the fit. The default is to output the log scaled predictions.

References

Miller, D.M. (1984) Reducing Transformation Bias in Curve Fitting The American Statistician. 38(2):124-126

Naihua, D. (1983) Smearing Estimate: A Nonparametric Retransformation Method Journal of the American Statistical Association, 78(383):605–610.

Examples

library(nlme)
  data(ralu.model)

back.transform <- function(y) exp(y + 0.5 * stats::var(y, na.rm=TRUE))
rmse = function(p, o){ sqrt(mean((p - o)^2)) } 

x = c("DEPTH_F", "HLI_F", "CTI_F", "cti", "ffp")
 
sidx <- sample(1:nrow(ralu.model), 100) 
  train <- ralu.model[sidx,]
  test <- ralu.model[-sidx,]
 
 # Specify constrained gravity model	
 ( gm <- gravity(y = "DPS", x = x, d = "DISTANCE", group = "FROM_SITE", 
                 data = train, ln = FALSE) )
  
( p <- predict(gm, test[,c(x, "DISTANCE")]) )
  rmse(back.transform(p), back.transform(ralu.model[,"DPS"][-sidx]))

# WIth model sigma-based back transformation
( p <- predict(gm, test[,c(x, "DISTANCE")], back.transform = "simple") )
( p <- predict(gm, test[,c(x, "DISTANCE")], back.transform = "Miller") )
( p <- predict(gm, test[,c(x, "DISTANCE")], back.transform = "Naihua") )

# Using grouped data
test <- nlme::groupedData(stats::as.formula(paste(paste("DPS", 1, sep = " ~ "), 
          "FROM_SITE", sep = " | ")), 
		  data = test[,c("DPS", "FROM_SITE", x, "DISTANCE")])

( p <- predict(gm, test, groups = "FROM_SITE") )
( y.hat <- back.transform(ralu.model[,"DPS"][-sidx]) )
    na.idx <- which(is.na(p))
  rmse(back.transform(p)[-na.idx], y.hat[-na.idx])

# Specify unconstrained gravity model (generally, not recommended)	
( gm <- gravity(y = "DPS", x = x, d = "DISTANCE", group = "FROM_SITE", 
                data = train, ln = FALSE, constrained=TRUE) )

( p <- predict(gm, test[,c(x, "DISTANCE")]) )
  rmse(back.transform(p), back.transform(ralu.model[,"DPS"][-sidx]))