Skip to content

Commit

Permalink
Regression bug fixes (#350)
Browse files Browse the repository at this point in the history 
* nTrials should be non-zero

* only estimates when there are predictors in best model

* update mahalanobis

* include try on confint call

* better error handling GLM log CI
  • Loading branch information
@JohnnyDoorn
JohnnyDoorn authored Nov 19, 2024
1 parent e0ce706 commit 5acdef4
Show file tree
Hide file tree
Showing 3 changed files with 41 additions and 24 deletions.
9 changes: 8 additions & 1 deletion R/commonglm.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -867,7 +867,14 @@

modelMatrix <- as.data.frame(model.matrix(model))
modelMatrix <- modelMatrix[colnames(modelMatrix) != "(Intercept)"]
influenceResData[["mahalanobis"]] <- mahalanobis(modelMatrix, center = colMeans(modelMatrix), cov = cov(modelMatrix))

if (ncol(modelMatrix) > 0) {
influenceResData[["mahalanobis"]] <- mahalanobis(modelMatrix, center = colMeans(modelMatrix), cov = cov(modelMatrix))
} else if (options[["mahalanobis"]]) {
influenceTable$addFootnote(
gettext("Mahalanobis distance cannot be computed for the intercept only model.")
)
}

if (options$residualCasewiseDiagnosticType == "cooksDistance")
index <- which(abs(influenceResData[["cooksDistance"]]) > options$residualCasewiseDiagnosticCooksDistanceThreshold)
Expand Down
21 changes: 14 additions & 7 deletions R/generalizedlinearmodel.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -102,21 +102,24 @@ GeneralizedLinearModelInternal <- function(jaspResults, dataset = NULL, options,
limits.min = 0,
limits.max = Inf,
exitAnalysisIfErrors = TRUE)

if (length(options$factors) != 0)
.hasErrors(dataset,
type = "factorLevels",
factorLevels.target = options$factors,
factorLevels.amount = '< 2',
exitAnalysisIfErrors = TRUE)

if (options[["family"]] == "bernoulli") {

if (length(levels(dataset[, options[["dependent"]]])) != 2)
.quitAnalysis(gettext("The Bernoulli family requires the dependent variable to be a factor with 2 levels."))

} else if (options[["family"]] == "binomial") {

if (any(dataset[, options[["weights"]]] == 0))
.quitAnalysis(gettext("The Binomial family requires the weights variable (i.e. total number of trials) to be non-zero."))

if (any(dataset[, options[["dependent"]]] < 0) || any(dataset[, options[["dependent"]]] > 1))
.quitAnalysis(gettext("The Binomial family requires the dependent variable (i.e. proportion of successes) to be between 0 and 1 (inclusive)."))

Expand Down Expand Up @@ -414,7 +417,11 @@ GeneralizedLinearModelInternal <- function(jaspResults, dataset = NULL, options,
rowNames <- rownames(modelSummary)

if (options[["coefficientCi"]]) {
coefCiSummary <- confint(fullModel, level = options[["coefficientCiLevel"]])
coefCiSummary <- try(confint(fullModel, level = options[["coefficientCiLevel"]]))
if (jaspBase::isTryError(coefCiSummary)) {
jaspResults[["estimatesTable"]]$setError("Confidence intervals not available for this model, try other predictors, families, or links.")
return()
}
if (length(rowNames) == 1) coefCiSummary <- matrix(coefCiSummary, ncol = 2)
} else {
coefCiSummary <- matrix(nrow = length(rowNames),
Expand Down Expand Up @@ -477,14 +484,14 @@ GeneralizedLinearModelInternal <- function(jaspResults, dataset = NULL, options,
.glmOutlierTable(jaspResults, dataset, options, ready, position = 8, residType = "standardized deviance")
.glmOutlierTable(jaspResults, dataset, options, ready, position = 8, residType = "studentized deviance")

.glmInfluenceTable(jaspResults[["diagnosticsContainer"]],
.glmInfluenceTable(jaspResults[["diagnosticsContainer"]],
jaspResults[["glmModels"]][["object"]][["fullModel"]],
dataset, options, ready, position = 9)
.regressionExportResiduals(jaspResults,

.regressionExportResiduals(jaspResults,
jaspResults[["glmModels"]][["object"]][["fullModel"]],
dataset, options, ready)

.glmMulticolliTable(jaspResults, dataset, options, ready, position = 10)

return()
Expand Down
35 changes: 19 additions & 16 deletions R/regressionlogistic.R
View file
Original file line number Diff line number Diff line change
Expand Up @@ -361,9 +361,9 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
rows <- vector("list", length(glmObj))

for (midx in seq_along(glmObj)) {

if (options$method == "enter")
.linregAddPredictorsInModelFootnote(jaspResults[["modelSummary"]],
.linregAddPredictorsInModelFootnote(jaspResults[["modelSummary"]],
options[["modelTerms"]][[midx]][["components"]], midx)
mObj <- glmObj[[midx]]
if (midx > 1) {
Expand Down Expand Up @@ -796,6 +796,9 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
mObj <- glmObj[[length(glmObj)]]
predictors <- character(0)

if (is.null(options$modelTerms[[length(glmObj)]][["components"]])) {
return()
}
mComponents <- options$modelTerms[[length(glmObj)]][["components"]]
predictors <- unlist(sapply(mComponents, function(x) if(length(x) == 1) x))

Expand Down Expand Up @@ -904,7 +907,7 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
container[["subContainer"]] <- createJaspContainer(gettext("Independent - predicted plots"))
subcontainer <- container[["subContainer"]]
if(!ready){
subcontainer[["placeholder"]] <- createJaspPlot(width = 480, height = 320, dependencies = c("independentVsPredictedPlot",
subcontainer[["placeholder"]] <- createJaspPlot(width = 480, height = 320, dependencies = c("independentVsPredictedPlot",
"independentVsPredictedPlotIncludeInteractions",
"independentVsPredictedPlotUseLogit"))
return()
Expand All @@ -914,32 +917,32 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
glmObj <- jaspResults[["glmRes"]][["object"]]
mObj <- glmObj[[length(glmObj)]]
mComponents <- options$modelTerms[[length(glmObj)]][["components"]]

if (options[["independentVsPredictedPlotIncludeInteractions"]]) {
predictors <- sapply(mComponents, function(x) if(length(x) < 3) x)
} else {
predictors <- sapply(mComponents, function(x) if(length(x) == 1) x)
}

predictions <- predict(mObj, type = "response")
if(options[["independentVsPredictedPlotUseLogit"]]) {
predictions <- log(predictions / (1 - predictions))
predictions <- log(predictions / (1 - predictions))
yName <- "Logit Predicted Probability"
} else {
yName <- "Predicted Probability"
}

for (pred in predictors) {

facPredictorIndex <- which(pred %in% options[["factors"]])

for (i in seq_along(pred)) {

predictorLogitPlot <- createJaspPlot(title = paste(c(pred[i], pred[-i]), collapse = " \u273B "), width = 480, height = 320)
predictorLogitPlot$dependOn(c("independentVsPredictedPlot",
predictorLogitPlot$dependOn(c("independentVsPredictedPlot",
"independentVsPredictedPlotIncludeInteractions",
"independentVsPredictedPlotUseLogit"))

binContVar <- FALSE
if (length(pred) == 1) {
groupVar <- groupName <- NULL
Expand Down Expand Up @@ -967,13 +970,13 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
yName = yName,
addSmooth = TRUE,
addSmoothCI = TRUE,
plotAbove = "none",
plotAbove = "none",
plotRight = "none",
showLegend = length(pred) > 1,
legendTitle = if (binContVar) paste0(groupName,"_binned") else groupName,
smoothCIValue = 0.95,
forceLinearSmooth = options[["independentVsPredictedPlotUseLogit"]]))

if(isTryError(p))
predictorLogitPlot$setError(.extractErrorMessage(p))
else {
Expand All @@ -982,10 +985,10 @@ RegressionLogisticInternal <- function(jaspResults, dataset = NULL, options, ...
predictorLogitPlot$plotObject <- p
}
subcontainer[[paste0(indepVar, groupName)]] <- predictorLogitPlot

}
}

return()
}

Expand Down

0 comments on commit 5acdef4

Please sign in to comment.