diff --git a/R/MixedModelsCommon.R b/R/MixedModelsCommon.R
index c8bc4f77..aac3df89 100644
--- a/R/MixedModelsCommon.R
+++ b/R/MixedModelsCommon.R
@@ -411,14 +411,92 @@
 
   return(dataset)
 }
-.mixedInterceptML   <- function(formula, dataset, type, family = NULL) {
+.mmCreateOptimizerControl <- function(type, options) {
+  # Create optimizer control objects based on user settings
+
+  # Build control arguments from user options
+  control_args <- list()
+
+  # Set optimizer method
+  if (options$optimizerMethod != "default") {
+    control_args$optimizer <- options$optimizerMethod
+  }
+
+  # Build optimizer-specific control list
+  optCtrl <- list()
+
+  if (options$optimizerMethod == "Nelder_Mead") {
+    # Nelder-Mead specific options
+    if (!is.null(options$nelderMeadMaxfun) && options$nelderMeadMaxfun > 0) {
+      optCtrl$maxfun <- options$nelderMeadMaxfun
+    }
+    if (!is.null(options$nelderMeadFtolAbs) && options$nelderMeadFtolAbs > 0) {
+      optCtrl$FtolAbs <- options$nelderMeadFtolAbs
+    }
+    if (!is.null(options$nelderMeadFtolRel) && options$nelderMeadFtolRel > 0) {
+      optCtrl$FtolRel <- options$nelderMeadFtolRel
+    }
+    if (!is.null(options$nelderMeadXtolRel) && options$nelderMeadXtolRel > 0) {
+      optCtrl$XtolRel <- options$nelderMeadXtolRel
+    }
+  } else if (options$optimizerMethod == "bobyqa") {
+    # bobyqa specific options
+    if (!is.null(options$bobyqaNpt) && options$bobyqaNpt > 0) {
+      optCtrl$npt <- options$bobyqaNpt
+    }
+    if (!is.null(options$bobyqaRhobeg) && options$bobyqaRhobeg > 0) {
+      optCtrl$rhobeg <- options$bobyqaRhobeg
+    }
+    if (!is.null(options$bobyqaRhoend) && options$bobyqaRhoend > 0) {
+      optCtrl$rhoend <- options$bobyqaRhoend
+    }
+    if (!is.null(options$bobyqaMaxfun) && options$bobyqaMaxfun > 0) {
+      optCtrl$maxfun <- options$bobyqaMaxfun
+    }
+  } else if (options$optimizerMethod == "nlminb") {
+    # nlminb specific options
+    if (!is.null(options$nlminbTol) && options$nlminbTol > 0) {
+      optCtrl$tol <- options$nlminbTol
+    }
+    if (!is.null(options$nlminbRelTol) && options$nlminbRelTol > 0) {
+      optCtrl$relTol <- options$nlminbRelTol
+    }
+    if (!is.null(options$nlminbMaxit) && options$nlminbMaxit > 0) {
+      optCtrl$maxit <- options$nlminbMaxit
+    }
+  } else {
+    # Default and BFGS: use generic options
+    if (!is.null(options$optimizerMaxIter) && options$optimizerMaxIter > 0) {
+      optCtrl$maxit <- options$optimizerMaxIter
+    }
+    if (!is.null(options$optimizerTolerance) && options$optimizerTolerance > 0) {
+      optCtrl$reltol <- options$optimizerTolerance
+    }
+  }
+
+  if (length(optCtrl) > 0) {
+    control_args$optCtrl <- optCtrl
+  }
+
+  # Create appropriate control object
+  if (type == "LMM") {
+    return(do.call(lme4::lmerControl, control_args))
+  } else if (type == "GLMM") {
+    return(do.call(lme4::glmerControl, control_args))
+  }
+}
+
+.mixedInterceptML   <- function(formula, dataset, type, family = NULL, options = NULL) {
   # this is a simple function to fit a mixed-effects model with a fixed intercept only
   # because afex does not allow those models for GLMMs (or LMMs with LRT/PB)
+  lmControl <<- .mmCreateOptimizerControl(type, options)
+
   if (type == "LMM") {
     fit <- lmerTest::lmer(
       formula         = formula,
       data            = dataset,
-      REML            = FALSE
+      REML            = FALSE,
+      control         = lmControl
     )
   } else if (type == "GLMM") {
     fit <- lme4::glmer(
@@ -474,7 +552,6 @@
   return(added)
 }
 .mmFitModel      <- function(jaspResults, dataset, options, type = "LMM") {
-
   if (!is.null(jaspResults[["mmModel"]]))
     return()
 
@@ -498,13 +575,19 @@
   # specify contrasts
   dataset <- .mmSetContrasts(dataset, options)
 
+  # the control arguments needs to be assigned outside of the call because
+  # forwarding the call crashes afex
+  lmControl <<- .mmCreateOptimizerControl(type, options)
+
   if (type == "LMM") {
     if (.isInterceptML(options))
       model <- try(
         .mixedInterceptML(
           formula         = as.formula(modelFormula$modelFormula),
           data            = dataset,
-          type            = "LMM"
+          type            = "LMM",
+          options         = options,
+          control         = lmControl
         ))
     else
       model <- try(
@@ -515,7 +598,8 @@
           method          = .mmGetTestMethod(options),
           test_intercept  = .mmGetTestIntercept(options),
           args_test       = list(nsim = options$bootstrapSamples),
-          check_contrasts = FALSE
+          check_contrasts = FALSE,
+          control         = lmControl
         ))
   } else if (type == "GLMM") {
     # needs to be evaluated in the global environment
@@ -540,7 +624,8 @@
           args_test       = list(nsim = options$bootstrapSamples),
           check_contrasts = FALSE,
           family          = glmmFamily,
-          weights         = glmmWeight
+          weights         = glmmWeight,
+          control         = lmControl
         ))
     } else {
       if (.isInterceptML(options))
@@ -549,7 +634,9 @@
             formula         = as.formula(modelFormula$modelFormula),
             data            = dataset,
             family          = glmmFamily,
-            type            = "GLMM"
+            type            = "GLMM",
+            options         = options,
+            control         = lmControl
           ))
       else
         model <- try(
@@ -562,7 +649,8 @@
             args_test       = list(nsim = options$bootstrapSamples),
             check_contrasts = FALSE,
             #start           = start,
-            family          = glmmFamily
+            family          = glmmFamily,
+            control         = lmControl
         ))
     }
   }
diff --git a/R/MixedModelsMessages.R b/R/MixedModelsMessages.R
index 721809fb..8b7d5f8a 100644
--- a/R/MixedModelsMessages.R
+++ b/R/MixedModelsMessages.R
@@ -237,7 +237,7 @@
   else if (grepl("Downdated VtV is not positive definite", error))
     return(gettext("The optimizer failed to find a solution. Probably due to scaling issues quasi-separation in the data. Try rescaling or removing some of the predictors."))
   else if (grepl("did not converge in (maxit) iterations", error))
-    return(gettext("The optimizer failed to find a solution in the specified number of iterations. (JASP currently does not support modifying the optimizer settings.)"))
+    return(gettext("The optimizer failed to find a solution in the specified number of iterations. Try adjusting the optimizer settings in the Advanced Options section."))
   else if (grepl("unexpected symbol", error)) # triggered by right hand side formula larger than 500 characters -- the maximum length
     return(gettext("The model formula is probably too long. Try shortening variable names."))
   else
diff --git a/inst/qml/MixedModelsGLMM.qml b/inst/qml/MixedModelsGLMM.qml
index 98de4c3e..012bde9d 100644
--- a/inst/qml/MixedModelsGLMM.qml
+++ b/inst/qml/MixedModelsGLMM.qml
@@ -419,4 +419,6 @@ Form {
 		}
 	}
 
+	MM.MixedModelsAdvanced {}
+
 }
diff --git a/inst/qml/MixedModelsLMM.qml b/inst/qml/MixedModelsLMM.qml
index f081e783..74ddcd1b 100755
--- a/inst/qml/MixedModelsLMM.qml
+++ b/inst/qml/MixedModelsLMM.qml
@@ -317,4 +317,6 @@ Form {
 		}
 	}
 
+	MM.MixedModelsAdvanced {}
+
 }
diff --git a/inst/qml/common/MixedModelsAdvanced.qml b/inst/qml/common/MixedModelsAdvanced.qml
new file mode 100644
index 00000000..bb34f646
--- /dev/null
+++ b/inst/qml/common/MixedModelsAdvanced.qml
@@ -0,0 +1,216 @@
+//
+// Copyright (C) 2013-2020 University of Amsterdam
+//
+// This program is free software: you can redistribute it and/or modify
+// it under the terms of the GNU Affero General Public License as
+// published by the Free Software Foundation, either version 3 of the
+// License, or (at your option) any later version.
+//
+// This program 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 Affero General Public License for more details.
+//
+// You should have received a copy of the GNU Affero General Public
+// License along with this program.  If not, see
+// <http://www.gnu.org/licenses/>.
+//
+import QtQuick
+import QtQuick.Layouts
+import JASP
+import JASP.Controls
+
+Section
+{
+	title:		qsTr("Advanced")
+	expanded:	false
+
+	Group
+	{
+		title: qsTr("Optimizer Settings")
+
+		DropDown
+		{
+			id:		optimizerDropdown
+			name:	"optimizerMethod"
+			label:	qsTr("Optimizer")
+			info:	qsTr("Optimization algorithm used for parameter estimation")
+			values:
+			[
+				{ label: qsTr("Default"),		value: "default"},
+				{ label: "nlminb",				value: "nlminb"},
+				{ label: "BFGS",				value: "BFGS"},
+				{ label: "Nelder-Mead",			value: "Nelder_Mead"},
+				{ label: "bobyqa",				value: "bobyqa"}
+			]
+		}
+
+		// Nelder-Mead specific options
+		IntegerField
+		{
+			name:			"nelderMeadMaxfun"
+			label:			qsTr("Maximum function evaluations")
+			info:			qsTr("Maximum number of function evaluations allowed")
+			defaultValue:	10000
+			min:			100
+			max:			1000000
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "Nelder_Mead"
+		}
+
+		DoubleField
+		{
+			name:			"nelderMeadFtolAbs"
+			label:			qsTr("Absolute function tolerance")
+			info:			qsTr("Absolute tolerance on change in function values")
+			defaultValue:	1e-5
+			min:			1e-12
+			max:			1e-1
+			decimals:		12
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "Nelder_Mead"
+		}
+
+		DoubleField
+		{
+			name:			"nelderMeadFtolRel"
+			label:			qsTr("Relative function tolerance")
+			info:			qsTr("Relative tolerance on change in function values")
+			defaultValue:	1e-15
+			min:			1e-20
+			max:			1e-5
+			decimals:		20
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "Nelder_Mead"
+		}
+
+		DoubleField
+		{
+			name:			"nelderMeadXtolRel"
+			label:			qsTr("Relative parameter tolerance")
+			info:			qsTr("Relative tolerance on change in parameter values")
+			defaultValue:	1e-7
+			min:			1e-15
+			max:			1e-3
+			decimals:		15
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "Nelder_Mead"
+		}
+
+		// bobyqa specific options
+		IntegerField
+		{
+			name:			"bobyqaNpt"
+			label:			qsTr("Number of interpolation points")
+			info:			qsTr("Number of points used to approximate the objective function")
+			defaultValue:	0
+			min:			0
+			max:			100
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "bobyqa"
+		}
+
+		DoubleField
+		{
+			name:			"bobyqaRhobeg"
+			label:			qsTr("Initial trust region radius")
+			info:			qsTr("Initial value of the trust region radius")
+			defaultValue:	0
+			min:			0
+			max:			10
+			decimals:		6
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "bobyqa"
+		}
+
+		DoubleField
+		{
+			name:			"bobyqaRhoend"
+			label:			qsTr("Final trust region radius")
+			info:			qsTr("Final value of the trust region radius")
+			defaultValue:	0
+			min:			0
+			max:			1
+			decimals:		10
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "bobyqa"
+		}
+
+		IntegerField
+		{
+			name:			"bobyqaMaxfun"
+			label:			qsTr("Maximum function evaluations")
+			info:			qsTr("Maximum number of function evaluations allowed")
+			defaultValue:	10000
+			min:			100
+			max:			1000000
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "bobyqa"
+		}
+
+		// nlminb specific options
+		DoubleField
+		{
+			name:			"nlminbTol"
+			label:			qsTr("Tolerance")
+			info:			qsTr("Convergence tolerance")
+			defaultValue:	1e-6
+			min:			1e-12
+			max:			1e-3
+			decimals:		12
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "nlminb"
+		}
+
+		DoubleField
+		{
+			name:			"nlminbRelTol"
+			label:			qsTr("Relative tolerance")
+			info:			qsTr("Relative convergence tolerance")
+			defaultValue:	1e-10
+			min:			1e-20
+			max:			1e-5
+			decimals:		20
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "nlminb"
+		}
+
+		IntegerField
+		{
+			name:			"nlminbMaxit"
+			label:			qsTr("Maximum iterations")
+			info:			qsTr("Maximum number of iterations")
+			defaultValue:	10000
+			min:			100
+			max:			1000000
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "nlminb"
+		}
+
+		// Default and BFGS options (generic fallback)
+		IntegerField
+		{
+			name:			"optimizerMaxIter"
+			label:			qsTr("Maximum iterations")
+			info:			qsTr("Maximum number of iterations for the optimizer")
+			defaultValue:	10000
+			min:			100
+			max:			1000000
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "default" || optimizerDropdown.currentValue === "BFGS"
+		}
+
+		DoubleField
+		{
+			name:			"optimizerTolerance"
+			label:			qsTr("Convergence tolerance")
+			info:			qsTr("Convergence tolerance for parameter estimates")
+			defaultValue:	1e-6
+			min:			1e-12
+			max:			1e-3
+			decimals:		12
+			fieldWidth:		80 * jaspTheme.uiScale
+			visible:		optimizerDropdown.currentValue === "default" || optimizerDropdown.currentValue === "BFGS"
+		}
+	}
+}
\ No newline at end of file
diff --git a/tests/testthat/test-mixedmodelsglmm.R b/tests/testthat/test-mixedmodelsglmm.R
index cc160fdb..e2e2109d 100644
--- a/tests/testthat/test-mixedmodelsglmm.R
+++ b/tests/testthat/test-mixedmodelsglmm.R
@@ -973,3 +973,38 @@ context("Generalized Linear Mixed Models")
     jaspTools::expect_equal_plots(testPlot, "plot-glmm-5")
   })
 }
+
+### Test optimizer options for GLMM
+{
+  test_that("GLMM optimizer options can be set without errors", {
+    options <- jaspTools::analysisOptions("MixedModelsGLMM")
+    options$dependent <- "dependent"
+    options$fixedEffects <- list(list(components = "factor1"))
+    options$randomEffects <- list(list(
+      randomComponents = list(list(randomSlopes = FALSE, value = "grouping"))
+    ))
+    options$family <- "binomial"
+    options$link <- "logit"
+    
+    # Test custom optimizer settings
+    options$optimizerMethod <- "bobyqa"
+    options$optimizerMaxIter <- 5000
+    options$optimizerMaxFunEvals <- 50000
+    options$optimizerTolerance <- 1e-8
+    options$optimizerCheckConv <- TRUE
+    
+    # Create simple test dataset for binomial GLMM
+    dataset <- data.frame(
+      dependent = rbinom(40, 1, 0.5),
+      factor1 = factor(rep(c("A", "B"), each = 20)),
+      grouping = factor(rep(1:4, each = 10))
+    )
+    
+    # Should not error during options processing
+    results <- jaspTools::runAnalysis("MixedModelsGLMM", dataset = dataset, options)
+    
+    # Basic check that analysis ran and produced some output
+    expect_true(!is.null(results))
+    expect_true(length(results$results) > 0)
+  })
+}
diff --git a/tests/testthat/test-mixedmodelslmm.R b/tests/testthat/test-mixedmodelslmm.R
index 81edf48d..9975ff13 100644
--- a/tests/testthat/test-mixedmodelslmm.R
+++ b/tests/testthat/test-mixedmodelslmm.R
@@ -1067,3 +1067,35 @@ context("Linear Mixed Models")
   })
 }
 
+### Test optimizer options
+{
+  test_that("Optimizer options can be set without errors", {
+    options <- jaspTools::analysisOptions("MixedModelsLMM")
+    options$dependent <- "Variable4"
+    options$fixedEffects <- list(list(components = "Variable1"))
+    options$randomEffects <- list(list(
+      randomComponents = list(list(randomSlopes = FALSE, value = "Variable2"))
+    ))
+    
+    # Test custom optimizer settings
+    options$optimizerMethod <- "nlminb"
+    options$optimizerMaxIter <- 5000
+    options$optimizerMaxFunEvals <- 50000
+    options$optimizerTolerance <- 1e-8
+    options$optimizerCheckConv <- TRUE
+    
+    dataset <- data.frame(
+      Variable1 = factor(rep(c("A", "B"), each = 20)),
+      Variable2 = factor(rep(1:4, each = 10)),
+      Variable4 = rnorm(40)
+    )
+    
+    # Should not error during options processing
+    results <- jaspTools::runAnalysis("MixedModelsLMM", dataset = dataset, options)
+    
+    # Basic check that analysis ran and produced some output
+    expect_true(!is.null(results))
+    expect_true(length(results$results) > 0)
+  })
+}
+