diff --git a/Project.toml b/Project.toml
index 09c9c230..15d546f0 100644
--- a/Project.toml
+++ b/Project.toml
@@ -1,9 +1,10 @@
 name = "ControlSystemIdentification"
 uuid = "3abffc1c-5106-53b7-b354-a47bfc086282"
 authors = ["baggepinnen <baggepinnen@gmail.com>"]
-version = "2.11.8"
+version = "2.11.9"
 
 [deps]
+ADTypes = "47edcb42-4c32-4615-8424-f2b9edc5f35b"
 ChainRules = "082447d4-558c-5d27-93f4-14fc19e9eca2"
 ComponentArrays = "b0b7db55-cfe3-40fc-9ded-d10e2dbeff66"
 ControlSystemsBase = "aaaaaaaa-a6ca-5380-bf3e-84a91bcd477e"
@@ -35,6 +36,7 @@ LeastSquaresOptim = "0fc2ff8b-aaa3-5acd-a817-1944a5e08891"
 ControlSystemIdentificationLSOptExt = "LeastSquaresOptim"
 
 [compat]
+ADTypes = "1"
 ChainRules = "1.58"
 ComponentArrays = "0.14, 0.15"
 ControlSystemsBase = "1.8"
@@ -50,7 +52,7 @@ LinearAlgebra = "1.6"
 LowLevelParticleFilters = "2.0, 3.0"
 MatrixEquations = "1.1, 2.0"
 MonteCarloMeasurements = "1.0"
-Optim = "<0.20, 0.20, 0.21, 0.22, 1.0"
+Optim = "2"
 Parameters = "<0.13"
 QuadGK = "2.4"
 Random = "1.6"
diff --git a/docs/Project.toml b/docs/Project.toml
index 92a840b0..c500efb3 100644
--- a/docs/Project.toml
+++ b/docs/Project.toml
@@ -16,7 +16,7 @@ SeeToDee = "1c904df7-48cd-41e7-921b-d889ed9a470c"
 StaticArrays = "90137ffa-7385-5640-81b9-e52037218182"
 
 [compat]
-DataInterpolations = "7"
+DataInterpolations = "8"
 Documenter = "1"
 LowLevelParticleFilters = "3.3.5"
 PrettyTables = "2"
diff --git a/src/ControlSystemIdentification.jl b/src/ControlSystemIdentification.jl
index 88282159..a13b0a9b 100644
--- a/src/ControlSystemIdentification.jl
+++ b/src/ControlSystemIdentification.jl
@@ -1,6 +1,7 @@
 module ControlSystemIdentification
 
-using ComponentArrays,
+using ADTypes,
+    ComponentArrays,
     ControlSystemsBase,
     DelimitedFiles,
     FFTW,
diff --git a/src/arx.jl b/src/arx.jl
index e671290e..d492a107 100644
--- a/src/arx.jl
+++ b/src/arx.jl
@@ -739,9 +739,9 @@ function tfest(
         end + 0.001abs2(a[1]-1) # Promote a monic polynomial for regularization (could enforce it without loss of generality.
     end
 
-    res = Optim.optimize(loss, ComponentVector(p0), opt, opts, autodiff = :forward)
+    res = Optim.optimize(loss, ComponentVector(p0), opt, opts; autodiff = AutoForwardDiff())
     if refine
-        res = Optim.optimize(loss2, res.minimizer, opt, opts, autodiff = :forward)
+        res = Optim.optimize(loss2, res.minimizer, opt, opts; autodiff = AutoForwardDiff())
     end
     a0 = res.minimizer.a[1]
     tf(res.minimizer.b ./ a0, res.minimizer.a ./ a0)
diff --git a/src/pem.jl b/src/pem.jl
index 781c867d..bbffd6a9 100644
--- a/src/pem.jl
+++ b/src/pem.jl
@@ -53,7 +53,7 @@ using Optim, Optim.LineSearches
         optimizer = BFGS(
             linesearch = LineSearches.BackTracking(),
         ),
-        autodiff = :forward,
+        autodiff = AutoForwardDiff(),
         store_trace = true,
         show_trace  = true,
         show_every  = 50,
@@ -90,7 +90,7 @@ If a manually provided initial guess `sys0`, this must also be scaled appropriat
 - `focus`: `prediction` or `:simulation`. If `:simulation`, the `K` matrix will be zero.
 - `h`: Prediction horizon for the prediction error filter. Large values of `h` makes the problem computationally expensive. As `h` approaches infinity, the problem approaches the `focus = :simulation` case.
 - `optimizer`: One of Optim's optimizers
-- `autodiff`: Whether or not to use forward-mode AD to compute gradients. `:forward` (default) for forward-mode AD, or `:finite` for finite differences.
+- `autodiff`: Whether or not to use forward-mode AD to compute gradients. `AutoForwardDiff()` (default) for forward-mode AD, or `AutoFiniteDiff()` for finite differences.
 - `metric`: The metric used to measure residuals. Try, e.g., `abs` for better resistance to outliers.
 The rest of the arguments are related to `Optim.Options`.
 - `regularizer`: A function of the parameter vector and the corresponding `PredictionStateSpace/StateSpace` system that can be used to regularize the estimate.
@@ -175,7 +175,7 @@ function newpem(
         # alphaguess = LineSearches.InitialStatic(alpha = 0.95),
         linesearch = LineSearches.BackTracking(),
     ),
-    autodiff = :forward,
+    autodiff = AutoForwardDiff(),
     store_trace = true,
     show_trace = true,
     show_every = 50,
@@ -764,8 +764,8 @@ function structured_pem(
         optimizer,
         Optim.Options(;
             store_trace, show_trace, show_every, iterations, allow_f_increases,
-            time_limit, x_abstol, f_abstol, g_tol, f_calls_limit, g_calls_limit),
-        autodiff = :forward,
+            time_limit, x_abstol, f_abstol, g_tol, f_calls_limit, g_calls_limit);
+        autodiff = AutoForwardDiff(),
     )
     sys_opt, K_opt, x0_opt = inner_constructor(res.minimizer, constructor, d.Ts)
     sysp_opt = PredictionStateSpace(