R/stat.glmm.no.first.model.R

Defines functions stat.glmm.no.first.model

Documented in stat.glmm.no.first.model

# Copyright (C) 2018  Sebastian Sosa, Ivan Puga-Gonzalez, Hu Feng He,Peng Zhang, Xiaohua Xie, Cédric Sueur
#
# This file is part of Animal Network Toolkit Software (ANTs).
#
# ANT is a free software: you can redistribute it and/or modify it
# under the terms of the GNU General Public License as published by
# the Free Software Foundation, either version 2 of the License, or
# (at your option) any later version.
#
# ANT is distributed in the hope that it will be useful, but
# WITHOUT ANY WARRANTY; without even the implied warranty of
# MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
# GNU General Public License for more details.

#' @title Extracts statistical measures of interest in Generalized Linear Mixed Models
#' @description Performs Generalized Linear Mixed Models tests
#' @param ant an output of ANT function \code{perm.net.nl} with random factor stated, or output of ANT 'met' categories functions in which output of ANT functions \code{perm.ds.focal}, \code{perm.ds.grp} or \code{perm.net.lk} where multiple matrices have been used.
#' @param formula two-sided linear formula object describing both the fixed-effects and random-effects part of the model, with the response on the left of a ~ operator and the terms, separated by + operators, on the right. Random-effects terms are distinguished by vertical bars (|) separating expressions for design matrices from grouping factors. Two vertical bars (||) can be used to specify multiple uncorrelated random effects for the same grouping variable. (Because of the way it is implemented, the ||-syntax works only for design matrices containing numeric (continuous) predictors; to fit models with independent categorical effects, see dummy or the lmer_alt function from the afex package.).
#' @param family a GLM family, see \code{\link{glm}} and \code{\link{family}}.
#' @param oda the original data frame of associations when argument ant is obtained with perm.ds.grp or perm.ds.focal ANT functions.
#' @param progress a boolean indicating the visualization of the permutation process.
#' @param odf the original data frame
#' @param ... Extra arguments for \code{lmer} or \code{glmer} function only.
#' @details GLMM with permutation data.
#' @return Returns a list of 3 elements :
#' \itemize{
#' \item An object of class \code{\link{merMod}} (more specifically, an object of subclass lmerMod or glmerMod), for which many methods are available (e.g. methods(class="merMod")).
#' \item A data frame if the estimates of the permuted models.
#' \item A vector of integers indicating the permutations that returned model errors or warnings (e.g. model convergence issues) and for which new permutations were done.
#' }
#' @seealso \code{\link{lmer}} or \code{\link{glmer}}
#' @author Sebastian Sosa, Ivan Puga-Gonzalez.

stat.glmm.no.first.model <- function(ant, formula, family, oda = NULL, progress = TRUE, odf, ...) {
  tmp <- NULL
  if (is.character(family)) {
    fam <- family
  }
  else {
    if (attributes(family)$class == "family") {
      family <- family
      fam <- family$family
    }
    else {
      stop("Argument family is not a character or a family function.")
    }
  }
  
  if (attributes(ant)$ANT == "ANT data stream group sampling multiple matrices") {
    if (is.null(oda)) {
      stop("Argument oda cannot be NULL when argument ant is obtained with perm.ds.grp or perm.ds.focal ANT functions")
    }
    
    if (fam == "gaussian") {
      # Parametrisation in case of error or warning -------------------------------------------------
      # finding node metrics to permute in case of warning or error along GLMM on permuted data
      arguments <- all.vars(formula)
      metrics <- c(
        "degree", "outdegree", "indegree", "strength", "outstrength", "instrength", "affinityB", "affinity", "affinityW",
        "disparity", "indisparity", "outdisparity", "eigenB", "eigenU", "outeigen", "ineigen", "eigenW", "eigen", "lpB", 
        "lpW", "reach", "riB", "riW", "ri"
      )
      
      target.metrics <- metrics[metrics %in% arguments]
      
      # Removing node metrics from original data frame
      odf <- odf[, -c(df.col.findId(odf, target.metrics))]
      
      # Finding scan and control factor to redo data stream permutation
      Scan <- attributes(ant)$scan
      ctrlf <- attributes(ant)$ctrlf
      method <- attributes(ant)$method
      
      # GLMM along permutations ------------------------------------------
      tmp.env <- new.env()
      tmp.env$new.perm <- 0
      tmp.env$gbi <- NULL
      tmp.env$error <- NULL
      
      if (progress) {
        permuted <- lapply(seq_along(ant), function(i, ant, formula, odf, oda, target.metrics, Scan, ctrlf, method, fam, ...) {
          cat("  Processing permutation : ", i, "\r")
          attr(oda, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lmer(formula = formula, data = ant[[i]], ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(ant[[i]])$permutation
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) gbi or controlGBI 3) glmm estimates
            r <- redo.ds.grp(family = "gaussian", new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, oda = oda, odf = odf, target.metrics = target.metrics, formula = formula, Scan = Scan, method = method, ctrlf = ctrlf, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, ant = ant, formula, odf, oda, target.metrics, Scan = Scan, ctrlf = ctrlf, method = method, fam = fam, ...)
        cat("\n")
      }
      else {
        permuted <- lapply(seq_along(ant), function(i, ant, formula, odf, oda, target.metrics, Scan, ctrlf, method, fam, ...) {
          attr(oda, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = ant[[i]], ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(ant[[i]])$permutation
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) gbi or controlGBI 3) glmm estimates
            r <- redo.ds.grp(family = "gaussian", new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, oda = oda, odf = odf, target.metrics = target.metrics, formula = formula, Scan = Scan, method = method, ctrlf = ctrlf, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, ant = ant, formula, odf, oda, target.metrics, Scan = Scan, ctrlf = ctrlf, method = method, fam = fam, ...)
      }
    }
    if (fam != "gaussian") {
      # Parametrisation in case of error or warning -------------------------------------------------
      # finding node metrics to permute in case of warning or error along GLMM on permuted data
      arguments <- all.vars(formula)
      metrics <- c(
        "degree", "outdegree", "indegree", "strength", "outstrength", "instrength", "affinityB", "affinity", "affinityW",
        "disparity", "indisparity", "outdisparity", "eigenB", "eigenU", "outeigen", "ineigen", "eigenW", "eigen", "lpB", 
        "lpW", "reach", "riB", "riW", "ri"
      )
      
      target.metrics <- metrics[metrics %in% arguments]
      
      # Removing node metrics from original data frame
      odf <- odf[, -c(df.col.findId(odf, target.metrics))]
      
      # Finding scan and control factor to redo data stream permutation
      Scan <- attributes(ant)$scan
      ctrlf <- attributes(ant)$ctrlf
      method <- attributes(ant)$method
      
      
      # GLMM along permutations ------------------------------------------
      tmp.env <- new.env()
      tmp.env$new.perm <- 0
      tmp.env$gbi <- NULL
      tmp.env$error <- NULL
      if (progress) {
        permuted <- lapply(seq_along(ant), function(i, ant, formula, family, odf, oda, target.metrics, Scan, ctrlf, method, fam) {
          cat("  Processing permutation : ", i, "\r")
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = ant[[i]], family = family, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(ant[[i]])$permutation
            
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) gbi or controlGBI 3) glmm estimates
            r <- redo.ds.grp(family = family, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, oda = oda, odf = odf, target.metrics = target.metrics, formula = formula, Scan = Scan, method = method, ctrlf = ctrlf, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          result <- summary(r)$coefficients[, 1]
          return(result)
        }, ant = ant, formula, family, odf, oda, target.metrics, Scan, ctrlf, method, fam, ...)
        cat("\n")
      }
      else {
        permuted <- lapply(seq_along(ant), function(i, ant, formula, family, odf, oda, target.metrics, Scan, ctrlf, method, fam, ...) {
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = ant[[i]], family = family, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(ant[[i]])$permutation
            
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) gbi or controlGBI 3) glmm estimates
            r <- redo.ds.grp(family = family, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, oda = oda, odf = odf, target.metrics = target.metrics, formula = formula, Scan = Scan, method = method, ctrlf = ctrlf, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          result <- summary(r)$coefficients[, 1]
          return(result)
        }, ant = ant, formula, family, odf, oda, target.metrics, Scan, ctrlf, method, fam, ...)
      }
    }
    permuted <- do.call("rbind", permuted)
    result <- list("permutations" = permuted, "errors" = tmp.env$error)
    attr(result, "class") <- "ant glmm parallel"
    attr(result, "family") <- paste(family)
    attr(result, "formula") <- format(formula)
    cat("\n")
    return(result)
  }
  
  if (attributes(ant)$ANT == "ANT data stream focal sampling multiple matrices") {
    if (is.null(oda)) {
      stop("Argument oda cannot be NULL when argument ant is obtained with perm.ds.grp or perm.ds.focal ANT functions")
    }
    
    if (fam == "gaussian") {
      # Parametrisation in case of error or warning -------------------------------------------------
      # finding node metrics to permute in case of warning or error along GLMM on permuted data
      arguments <- all.vars(formula)
      metrics <- c(
        "degree", "outdegree", "indegree", "strength", "outstrength", "instrength", "affinityB", "affinity", "affinityW",
        "disparity", "indisparity", "outdisparity", "eigenB", "eigenU", "outeigen", "ineigen", "eigenW", "eigen", "lpB", 
        "lpW", "reach", "riB", "riW", "ri"
      )
      
      target.metrics <- metrics[metrics %in% arguments]
      
      # Removing node metrics from original data frame
      odf <- odf[, -c(df.col.findId(odf, target.metrics))]
      
      # Finding scan and control factor to redo data stream permutation
      focal <- attributes(ant)$focal
      ctrl <- attributes(ant)$ctrl
      alters <- attributes(ant)$alters
      index <- attributes(ant)$method
      
      # GLMM along permutations ------------------------------------------
      tmp.env <- new.env()
      tmp.env$new.perm <- 0
      tmp.env$gbi <- NULL
      tmp.env$gbi2 <- NULL
      tmp.env$error <- NULL
      
      if (progress == TRUE) {
        permuted <- lapply(ant, function(d, formula, odf, oda, target.metrics, focal, ctrl, alters, index, fam, ...) {
          cat("  Processing permutation : ", attributes(d)$permutation, "\r")
          attr(odf, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = d, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(d)$permutation
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) permuted data frame of associations 3) glmm estimates
            r <- redo.ds.focal.glmm(family = family, formula = formula, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, gbi2 = tmp.env$gbi2, oda = oda, odf = odf, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            tmp.env$gbi2 <- r[[3]]
            result <- r[[4]]
            return(result)
          }
          
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula = formula, odf = odf, oda = oda, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
        cat("\n")
      }
      else {
        permuted <- lapply(ant, function(d, formula, odf, oda, target.metrics, focal, ctrl, alters, index, fam, ...) {
          attr(odf, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = d, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(d)$permutation
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) permuted data frame of associations 3) glmm estimates
            r <- redo.ds.focal.glmm(family = family, formula = formula, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, gbi2 = tmp.env$gbi2, oda = oda, odf = odf, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
            tmp.env$new.perm <- r[[1]]
            tmp.env$error <- c(tmp.env$error, r[[1]])
            tmp.env$gbi <- r[[2]]
            tmp.env$gbi2 <- r[[3]]
            result <- r[[4]]
            return(result)
          }
          
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula = formula, odf = odf, oda = oda, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
      }
    }
    
    if (fam != "gaussian") {
      # Parametrisation in case of error or warning -------------------------------------------------
      # finding node metrics to permute in case of warning or error along GLMM on permuted data
      arguments <- all.vars(formula)
      metrics <- c(
        "degree", "outdegree", "indegree", "strength", "outstrength", "instrength", "affinityB", "affinity", "affinityW",
        "disparity", "indisparity", "outdisparity", "eigenB", "eigenU", "outeigen", "ineigen", "eigenW", "eigen", "lpB", 
        "lpW", "reach", "riB", "riW", "ri"
      )
      
      target.metrics <- metrics[metrics %in% arguments]
      
      # Removing node metrics from original data frame
      odf <- odf[, -c(df.col.findId(odf, target.metrics))]
      
      # Finding scan and control factor to redo data stream permutation
      focal <- attributes(ant)$focal
      ctrl <- attributes(ant)$ctrl
      alters <- attributes(ant)$alters
      
      
      # GLMM along permutations ------------------------------------------
      new.perm <- 0
      new.oda <- NULL
      er <- NULL
      if (progress == TRUE) {
        permuted <- lapply(ant, function(d, formula, family, progress = TRUE, odf, oda, target.metrics, focal, ctrl, alters, new.perm, new.oda, fam, ...) {
          cat("  Processing permutation : ", attributes(d)$permutation, "\r")
          attr(odf, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = d, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            er <- c(er, attributes(d)$permutation)
            cat("  Processing permutation ", attributes(d)$permutation, " resampling", "\r")
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(d)$permutation
            
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) permuted data frame of associations 3) glmm estimates
            r <- redo.ds.focal.glmm(family = family, formula = formula, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, gbi2 = tmp.env$gbi2, oda = oda, odf = odf, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
            
            new.perm <- r[[1]]
            new.oda <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula, family, progress = TRUE, odf = odf, oda = oda, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, new.perm = new.perm, new.oda = new.oda, fam = fam, ...)
        cat("\n")
      }
      else {
        permuted <- lapply(ant, function(d, formula, family, progress = TRUE, odf, oda, target.metrics, focal, ctrl, alters, new.perm, new.oda, fam, ...) {
          attr(odf, "permutation") <- 0
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = d, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            # redo a permutation on raw data
            er <- c(er, attributes(d)$permutation)
            cat("  Processing permutation ", attributes(d)$permutation, " resampling", "\r")
            # Giving to the original data frame of individual characteristics (odf) the permutation number where error or warning were found
            attr(odf, "permutation") <- attributes(d)$permutation
            
            # redo.ds.grp.first return 3 elements: 1) permutation index, 2) permuted data frame of associations 3) glmm estimates
            r <- redo.ds.focal.glmm(family = family, formula = formula, new.perm = tmp.env$new.perm, gbi = tmp.env$gbi, gbi2 = tmp.env$gbi2, oda = oda, odf = odf, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, index = index, fam = fam, ...)
            
            new.perm <- r[[1]]
            new.oda <- r[[2]]
            result <- r[[3]]
            return(result)
          }
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula, family, progress = TRUE, odf = odf, oda = oda, target.metrics = target.metrics, focal = focal, ctrl = ctrl, alters = alters, new.perm = new.perm, new.oda = new.oda, fam = fam, ...)
      }
    }
    
    permuted <- do.call("rbind", permuted)
    result <- list("permutations" = permuted, "errors" = tmp.env$error)
    attr(result, "class") <- "ant glmm parallel"
    attr(result, "family") <- paste(family)
    attr(result, "formula") <- format(formula)
    cat("\n")
    return(result)
  }
  
  if (attributes(ant)$ANT == "ANT node label permutation with random factors") {
    if (fam == "gaussian") {
      # Test along the list of permutation data ------------------------------------------------------------------------
      if (progress) {
        tmp.env <- new.env()
        tmp.env$error <- NULL
        ctrl <- attributes(ant)$rf
        labels <- attributes(ant)$labels
        
        permuted <- lapply(seq_along(ant), function(i, ant, formula, progress, ctrl, odf, labels, ...) {
          cat("  Processing permutation : ", attributes(ant[[i]])$permutation, "\r")
          
          r <- tryCatch(suppressWarnings(suppressMessages(lmer(formula = formula, data = ant[[i]], ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(tmp@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            tmp.env$error <- c(tmp.env$error, attributes(ant[[i]])$permutation)
            while (all(test) != TRUE) {
              newdf <- perm.redo(df = odf, labels = labels, ctrl = ctrl)
              r <- tryCatch(suppressWarnings(suppressMessages(lmer(formula = formula, data = newdf, ...))), error = identity)
              
              if (isS4(r)) {
                r2=with(r@optinfo$derivs,solve(Hessian,gradient))
                if(max(abs(r2))<0.001){test=TRUE}
                else{
                  test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
                }
              }
              if (is(r, "error")) {
                test <- FALSE
              }
              if (is(r, "warning")) {
                test <- FALSE
              }
            }
          }
          
          r <- summary(r)$coefficients[, 1]
          
          return(r)
        }, ant = ant, formula, progress = TRUE, ctrl = ctrl, odf = odf, labels = labels, ...)
        cat("\n")
      }
      else {
        tmp.env <- new.env()
        tmp.env$error <- NULL
        ctrl <- attributes(ant)$rf
        labels <- attributes(ant)$labels
        
        permuted <- lapply(seq_along(ant), function(i, ant, formula, progress, ctrl, odf, labels, ...) {
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = ant[[i]], ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            tmp.env$error <- c(tmp.env$error, attributes(ant[[i]])$permutation)
            while (all(test) != TRUE) {
              newdf <- perm.redo(df = odf, labels = labels, ctrl = ctrl)
              r <- tryCatch(suppressWarnings(suppressMessages(lme4::lmer(formula = formula, data = newdf, ...))), error = identity)
              
              if (isS4(r)) {
                r2=with(r@optinfo$derivs,solve(Hessian,gradient))
                if(max(abs(r2))<0.001){test=TRUE}
                else{
                  test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
                }
              }
              if (is(r, "error")) {
                test <- FALSE
              }
              if (is(r, "warning")) {
                test <- FALSE
              }
            }
          }
          
          r <- summary(r)$coefficients[, 1]
          
          return(r)
        }, ant = ant, formula, progress = TRUE, ctrl = ctrl, odf = odf, labels = labels, ...)
      }
    }
    
    if (fam != "gaussian") {
      # Test along the list of permutation data ------------------------------------------------------------------------
      if (progress == TRUE) {
        tmp.env <- new.env()
        tmp.env$error <- NULL
        
        permuted <- lapply(ant, function(d, formula, family, ctrl, labels, odf, ...) {
          cat("  Processing permutation : ", attributes(d)$permutation, "\r")
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = d, family = family, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            tmp.env$error <- c(tmp.env$error, attributes(d)$permutation)
          }
          
          while (all(test) != TRUE) {
            newdf <- perm.redo(df = odf, labels = labels, ctrl = ctrl)
            r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = newdf, family = family, ...))), error = identity)
            
            if (isS4(r)) {
              r2=with(r@optinfo$derivs,solve(Hessian,gradient))
              if(max(abs(r2))<0.001){test=TRUE}
              else{
                test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
              }
            }
            if (is(r, "error")) {
              test <- FALSE
            }
            if (is(r, "warning")) {
              test <- FALSE
            }
          }
          
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula = formula, family = family, ctrl = ctrl, labels, odf, ...)
        cat("\n")
      }
      else {
        tmp.env <- new.env()
        tmp.env$error <- NULL
        
        permuted <- lapply(ant, function(d, formula, family, ctrl, labels, odf, w, ...) {
          r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = d, family = family, ...))), error = identity)
          
          if (isS4(r)) {
            r2=with(r@optinfo$derivs,solve(Hessian,gradient))
            if(max(abs(r2))<0.001){test=TRUE}
            else{
              test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
            }
          }
          if (is(r, "error")) {
            test <- FALSE
          }
          if (is(r, "warning")) {
            test <- FALSE
          }
          
          if (all(test) != TRUE) {
            tmp.env$error <- c(tmp.env$error, attributes(d)$permutation)
          }
          
          while (all(test) != TRUE) {
            newdf <- perm.redo(df = odf, labels = labels, ctrl = ctrl)
            r <- tryCatch(suppressWarnings(suppressMessages(lme4::glmer(formula = formula, data = newdf, family = family, ...))), error = identity)
            if (isS4(r)) {
              r2=with(r@optinfo$derivs,solve(Hessian,gradient))
              if(max(abs(r2))<0.001){test=TRUE}
              else{
                test <- c(!is(r, "error"), !is(r, "warning"),r@optinfo$conv$opt == 0, length(r@optinfo$conv$lme4$messages) == 0, length(r@optinfo$warnings) == 0)
              }
            }
            if (is(r, "error")) {
              test <- FALSE
            }
            if (is(r, "warning")) {
              test <- FALSE
            }
          }
          
          r <- summary(r)$coefficients[, 1]
          return(r)
        }, formula = formula, family = family, ctrl = ctrl, labels, odf, ...)
      }
    }
    
    permuted <- do.call("rbind", permuted)
    result <- list("permutations" = permuted, "errors" = tmp.env$error)
    attr(result, "class") <- "ant glmm parallel"
    attr(result, "family") <- paste(fam)
    attr(result, "formula") <- format(formula)
    cat("\n")
    return(result)
  }
}

Try the ANTs package in your browser

Any scripts or data that you put into this service are public.

ANTs documentation built on July 3, 2022, 1:05 a.m.