This site contains statistical tools for analyzing data from satellite dive recorders, as described in:

Simpkins, M. A., K. L. Laidre and P. J. Heagerty. 2005. Multivariate regression of satellite-linked dive recorder data: simultaneous analysis of all bins. Marine Mammal Science in press.

Specifically, this site includes statistical code and descriptions for two integrated Splus 6.1 functions that use the output from the Splus glm function and estimate an empirical variance-covariance matrix using generalized estimating equations (GEE). Similar analyses can also be conducted using GEE options of the SAS procedure GENMOD, or using the “cluster(id)” option with the Poisson regression methods in STATA.

To utilize the functions below, one must first create a generalized linear model (GLM) using the Splus function glm. Simpkins et al. (2005) assumed complete independence among data points for the GLM of harbor seal diving behavior. The empirical variance estimation conducted within GEE, using the functions shown below, provides a valid large-sample estimate of the variance-covariance matrix for the regression estimates even when an incorrect correlation model (e.g., complete independence) is specified.

After the GLM analysis is complete, the resulting glm object can be used as input into the infjack.glm function shown below. This function, in turn, calls the estfun.glm function and estimates an empirical variance-covariance matrix for the GLM, using GEE techniques. The infjack.glm and estfun.glm functions shown below can be copied and pasted into Splus at the command line to create the functions within an active Splus 6.1 workspace. Text files infjack.txt and estfun.txt are available for download via links at the bottom of this page.

Interested parties should refer to Simpkins et al. (2005) for clarification regarding the use of these functions for analyzing SDR data. Further questions should be addressed to the authors directly.

Function 1: estfun.glm

estfun.glm <- function(glm.obj) 
{ 
## 
# Create X matrix from glm object # 
## 
      if(is.matrix(glm.obj$x)) 
           xmat <- glm.obj$x 
      else { 
           mf <- model.frame(glm.obj) 
           xmat <- model.matrix(terms(glm.obj), mf) 
      } 
## 
# Calculate deviance from original glm regression for each data point # 
## 
      output <- residuals(glm.obj, "working") * glm.obj$weights * xmat 
## 
# Output this point-by-point deviance # 
## 
      output 
}

infjack.glm <- function(glm.obj,groups) 
{ 
## 
# Run estfun.glm to get point-by-point deviance from glm regression # 
## 
      umat <- estfun.glm(glm.obj) 
## 
# Estimate variance within each individual's data record (i.e. within # 
# each "group") by summing the point-by-point deviance from the glm # 
# regression for each individual (group) # 
## 
      usum <- rowsum(umat,groups,reorder=F) 
## 
# Use variance-covariance matrix from glm regression (modelv below) # 
# and variance within each individual's data record (usum above) to # 
# estimate empirical variance-covariance matrix (output below) # 
## 
      modelv <- summary(glm.obj)$cov.unscaled 
      output <- modelv%*%(t(usum)%*%usum)%*%modelv 
## 
# Output empirical variance-covariance matrix # 
## 
      output 
}