[R] Fwd: Error message: object of type 'closure' is not subsettable

[Simrit Rattan]

hey everyone :),
Subject: Re: Error message: object of type 'closure' is not subsettable
I am writing a package which should calculate the binary logistic
regression.
The function itself work perfectly, but if I want to load the function from
my package it gives me the above mentioned error. I have tried a lot of
things, googled and so on, but I can not figure out what to change so it
works. I also found out the mistake lies in the first function
(logisticRegression)

logisticRegression <- function(x,y, threshold = 1e-10, maxIter = 100)
{
  calcPi <- function(x,betaCoef)
  {
    betaCoef <- as.matrix(betaCoef)
    return(exp(x%*%betaCoef)/(1 + exp(x%*%betaCoef)))
  }
  #initial guess for beta (mostly we start with 0)
  betaCoef <- rep(0, ncol(x))

  #some initial value which is bigger than the threshold for the loop
  initialValue <- 10000

  #count of iteration to make sure its not an infinite loop
  iterCount <- 0

  #iteration process (loop)
  while(initialValue > threshold) #convergence test
  {
    #calculates probabilities using the current estimate of beta
    p <- as.vector(calcPi(x, betaCoef))

    #calculates W which is needed for the Hessian
    W <- diag(p*(1-p))

    #calculates the gradient
    gradient <- t(x)%*%(y-p)

    #calculates the hessian
    hessian <- -(t(x)%*%W)%*%x

    #calculates beta
    betaCoef <- betaCoef - solve(hessian) %*% gradient

    #how much did we change beta?
    initialValue <- sum((solve(hessian) %*% gradient)^2)

    #to see if we reached the max # of iteration
    iterCount <- iterCount + 1
    if(iterCount > maxIter) {
      stop("This is not converging")

    }
  }
  #df
  Totaldf <- length(y) - 1
  Residualdf <- length(y) - length(betaCoef)
  degreesOfFreedom <- cbind(Totaldf, Residualdf)
  #fisher information, variance, standard error
  fisherInformation <- -hessian
  vCov <- solve(fisherInformation)
  stdError <- sqrt(diag(vCov))

  #wald test statistic
  waldTestStatistic <- betaCoef/stdError
  pValue <- 2 * pnorm(-abs(waldTestStatistic), mean = 0, sd = 1)

  #loglikelihood and deviance
  logLikelihoodFunction <-  y%*%x%*%betaCoef - sum(log(1+exp(x%*%betaCoef)))
  residualDeviance <- -2*logLikelihoodFunction

  #AIC
  akaikeIC <- residualDeviance + 2*ncol(betaCoef)
  #for R2
  cV <- 2/nrow(x) #for C&S R2

  #deviance residual
  sign1 <- replace(y, y == 0, -1)
  devianceResiduals <- sign1*sqrt(-2*(y*log(p) + (1-y)*log(1-p)))
  residualQuantiles <- summary(devianceResiduals)

  #fitted values
  fittedValues <- log(p/(1-p))
  #what to return
  coefficients <- betaCoef
  output <- list()
  output$degreesF <- Residualdf
  output$nulldf <- Totaldf
  output$fittedValues <- fittedValues
  output$devianceResiduals <- devianceResiduals
  output$residualQuantiles <- residualQuantiles
  output$coefficients <- coefficients
  output$vCov <- vCov
  output$coefMat <- cbind(coefficients,
                          stdError,
                          waldTestStatistic,
                          pValue)
  colnames(output$coefMat) <- c("Estimate",
                                "Std.Err",
                                "z value",
                                "Pr(>|z|)")
  output$logLikelihoodFunction <- logLikelihoodFunction
  output$residualDeviance <- residualDeviance
  output$akaikeIC <- akaikeIC
  output$iterCount <- iterCount -1
  output$p <- p
  output$cV <- cV
  output


}

logReg <- function(formula, data) {
  x <- model.matrix(formula, data)
  regressant <- model.frame(formula, data)
  y <- regressant[,1]
  constant <- matrix(rep(1, length(y)))
  nulllog <- logisticRegression(constant,y)
  nullDeviance <- nulllog$residualDeviance
  nullLogLikelihood <- nulllog$logLikelihoodFunction

  output <- list()
  output <- logisticRegression(x,y)
  output$formula <- formula
  output$nullDeviance <- nullDeviance
  output$nullLogLikelihood <- nullLogLikelihood
  output$call <- match.call()

  #define the S3 class

  class(output) <- "logReg"
  print(output)

}

logReg.default <- function(x, ...) UseMethod("logReg")

print.logReg <- function(x, ...) {
  cat("\nCall:\n")
  print(x$call)
  cat("\nCoefficients:\n")
  colnames(x$coefficients) <- ""
  print(t(x$coefficients), digits = 5)
  cat("\nDegrees of Freedom:",
  paste(x$nulldf, "Total (i.e. Null);", x$degreesF, "Residual"))
  cat("\nNull Deviance:", paste(round(x$nullDeviance, digits = 2)))
  cat("\nResidual Deviance:", paste(round(x$residualDeviance, digits = 2)))
  invisible(x)
}

Thanks for your help! :D




--
Sent from: http://r.789695.n4.nabble.com/R-help-f789696.html

	[[alternative HTML version deleted]]