diff --git a/src/api/estimators.jl b/src/api/estimators.jl
index baa284f..b3ac7d7 100644
--- a/src/api/estimators.jl
+++ b/src/api/estimators.jl
@@ -67,10 +67,10 @@ available_optlib(dre::DensityRatioEstimator) = available_optlib(typeof(dre))
 ## functions to be implemented by new estimators ##
 ###################################################
 
-_densratio(x_nu, x_de, dre::DensityRatioEstimator, optlib::Type{OptimizationLibrary}) = @error "not implemented"
+_densratio(x_nu, x_de, dre::DensityRatioEstimator, optlib::Type{<:OptimizationLibrary}) = _throw_not_implemented_error("densratio", dre, optlib)
 
-_densratiofunc(x_nu, x_de, dre::DensityRatioEstimator, optlib::Type{OptimizationLibrary}) = @error "not implemented"
+_densratiofunc(x_nu, x_de, dre::DensityRatioEstimator, optlib::Type{<:OptimizationLibrary}) = _throw_not_implemented_error("densratiofunc", dre, optlib)
 
-default_optlib(dre::Type{DensityRatioEstimator}) = @error "not implemented"
+default_optlib(dre::Type{<:DensityRatioEstimator}) = _throw_not_implemented_error("default_optlib", dre)
 
-available_optlib(dre::Type{DensityRatioEstimator}) = @error "not implemented"
+available_optlib(dre::Type{<:DensityRatioEstimator}) = _throw_not_implemented_error("available_optlib", dre)
diff --git a/src/api/fitters.jl b/src/api/fitters.jl
index e41dc25..12a6b53 100644
--- a/src/api/fitters.jl
+++ b/src/api/fitters.jl
@@ -31,5 +31,4 @@ fit(dre::Type{<:DensityRatioEstimator}, x_nu, x_de, fitter::EstimatorFitter; opt
 ## functions to be implemented by new fitters ##
 ################################################
 
-_fit(dre::Type{DensityRatioEstimator}, x_nu, x_de, fitter::EstimatorFitter, optlib::Type{OptimizationLibrary}) =
-  @error "not implemented"
+_fit(dre::Type{<:DensityRatioEstimator}, x_nu, x_de, fitter::EstimatorFitter, optlib::Type{<:OptimizationLibrary}) = _throw_not_fit_error(dre, fitter, optlib)
\ No newline at end of file
diff --git a/src/kliep.jl b/src/kliep.jl
index 02310bf..2c61b2f 100644
--- a/src/kliep.jl
+++ b/src/kliep.jl
@@ -60,4 +60,4 @@ end
 
 Return the coefficients of KLIEP basis expansion.
 """
-_kliep_coeffs(K_nu, K_de, dre::KLIEP, optlib::Type{<:OptimizationLibrary}) = @error "not implemented"
+_kliep_coeffs(K_nu, K_de, dre::KLIEP, optlib::Type{<:OptimizationLibrary}) = _throw_opt_error(dre, optlib)
diff --git a/src/kmm.jl b/src/kmm.jl
index bbd6c30..84da514 100644
--- a/src/kmm.jl
+++ b/src/kmm.jl
@@ -86,3 +86,5 @@ end
 default_optlib(dre::Type{<:KMM}) = JuMPLib
 
 available_optlib(dre::Type{<:KMM}) = [JuMPLib]
+
+_kmm_ratios(K, κ, dre::AbstractKMM, optlib::Type{<:OptimizationLibrary}) = _throw_opt_error(dre, optlib)
diff --git a/src/lsif.jl b/src/lsif.jl
index a98beb2..c23cac8 100644
--- a/src/lsif.jl
+++ b/src/lsif.jl
@@ -73,4 +73,4 @@ end
 
 Return the coefficients of LSIF basis expansion.
 """
-_lsif_coeffs(H, h, dre::LSIF, optlib::Type{<:OptimizationLibrary}) = @error "not implemented"
+_lsif_coeffs(H, h, dre::LSIF, optlib::Type{<:OptimizationLibrary}) = _throw_opt_error(dre, optlib)
diff --git a/src/utils.jl b/src/utils.jl
index cb82488..cac269f 100644
--- a/src/utils.jl
+++ b/src/utils.jl
@@ -61,3 +61,30 @@ It is compatible with
 - Zygote.jl (see lib/zygote.jl)
 """
 safe_diagm(mat, a) = a * I
+
+
+###################################################
+##   Functions and objects for throwing errors   ##
+###################################################
+
+OPTLIB_DICT = Dict(
+    "JuliaLib" => "Julia",
+    "OptimLib" => "Optim",
+    "ConvexLib" => "Convex",
+    "JuMPLib" => "JuMP"
+)
+
+function _throw_opt_error(dre::DensityRatioEstimator, optlib::Type{<:OptimizationLibrary})
+  dre_name = nameof(typeof(dre))
+  lib_name = OPTLIB_DICT[string(optlib)]
+  optlib_options = join(available_optlib(dre), ", ")
+  error("Attempted to compute $(dre_name) density ratios using (possibly default) optimization library $(optlib), but this library has either not been loaded or is not implemented for use with $(dre_name). Available options for `optlib`: $(optlib_options). If $(optlib) is contained within the available options, be sure to call `using $(lib_name)` before calling `densratio`, `densratiofunc`, or other functions fitting the $(dre_name) estimator.")
+end
+
+function _throw_not_implemented_error(func::String, dre::Type{<:DensityRatioEstimator})
+  dre_name = nameof(dre)
+  error("Attempted to call `$(func)($(dre_name))` but this function has not been implemented for density ratio estimator of type $(dre_name).")
+end
+
+_throw_not_fit_error(dre::Type{<:DensityRatioEstimator}, fitter::EstimatorFitter, optlib::Type{<:OptimizationLibrary}) = 
+  error("Attempted to `fit` estimator `$(nameof(dre))` using fitter `$(typeof(fitter))` with optimization library `optlib=$(optlib)`, but no `fit` function has been implemented for this combination of estimator, fitter, and optimization library.")
diff --git a/test/runtests.jl b/test/runtests.jl
index 04d8b33..affa385 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -2,13 +2,44 @@ using DensityRatioEstimation
 using Distributions
 using LinearAlgebra
 using Statistics
+using Test, StableRNGs
+using ReferenceTests, ImageIO
+import CairoMakie as Mke
+
+@testset "Throw errors: optlib not loaded" begin
+
+  @test_throws ErrorException densratio([0], [0], KLIEP())
+  @test_throws ErrorException densratiofunc([0], [0], KLIEP())
+  @test_throws ErrorException densratio([0], [0], LSIF())
+  @test_throws ErrorException densratiofunc([0], [0], LSIF())
+  @test_throws ErrorException densratio([0], [0], KMM(); optlib = JuMPLib)
+
+end
+
 using Optim
 using JuMP, Ipopt
 using Convex, ECOS
-using Test, StableRNGs
-using ReferenceTests, ImageIO
 
-import CairoMakie as Mke
+@testset "Throw errors: optlib not available" begin
+
+  @test_throws ErrorException densratio([0], [0], KLIEP(); optlib = JuMPLib)
+  @test_throws ErrorException densratio([0], [0], LSIF(); optlib = ConvexLib)
+  @test_throws ErrorException densratio([0], [0], KMM(); optlib = ConvexLib)
+
+end
+
+@testset "Throw errors: undefined functions" begin
+
+  @test_throws ErrorException fit(KMM, [0], [0], LCV((σ=[1.],b=[10])))
+
+  struct NewDRE <: DensityRatioEstimator end
+
+  @test_throws ErrorException densratio([0], [0], NewDRE())
+  @test_throws ErrorException densratiofunc([0], [0], NewDRE())
+  @test_throws ErrorException default_optlib(NewDRE)
+  @test_throws ErrorException available_optlib(NewDRE)
+
+end
 
 # environment settings
 isCI = "CI" ∈ keys(ENV)