diff --git a/docs/src/basics/AbstractSystem.md b/docs/src/basics/AbstractSystem.md
index e68a7bb94e..3aa43e6d25 100644
--- a/docs/src/basics/AbstractSystem.md
+++ b/docs/src/basics/AbstractSystem.md
@@ -69,6 +69,8 @@ Optionally, a system could have:
 Note that if you know a system is an `AbstractTimeDependentSystem` you could use `get_iv` to get the
 unique independent variable directly, rather than using `independent_variables(sys)[1]`, which is clunky and may cause problems if `sys` is an `AbstractMultivariateSystem` because there may be more than one independent variable. `AbstractTimeIndependentSystem`s do not have a method `get_iv`, and `independent_variables(sys)` will return a size-zero result for such. For an `AbstractMultivariateSystem`, `get_ivs` is equivalent.
 
+For the `parameters`, `unknowns`, `continuous_events`, and `discrete_events` accessors there are corresponding `parameters_toplevel`, `unknowns_toplevel`, `continuous_events_toplevel`, and `discrete_events_toplevel` accessors which work similarly, but ignore the content of subsystems. Furthermore, a `equations_toplevel` version of `equations` exists as well, however, it can only be applied to non-complete systems.
+
 A system could also have caches:
 
   - `get_jac(sys)`: The Jacobian of a system.
diff --git a/src/systems/abstractsystem.jl b/src/systems/abstractsystem.jl
index c0b3601e21..9e6e875283 100644
--- a/src/systems/abstractsystem.jl
+++ b/src/systems/abstractsystem.jl
@@ -1313,6 +1313,18 @@ function unknowns(sys::AbstractSystem)
     unique(nonunique_unknowns)
 end
 
+"""
+    unknowns_toplevel(sys::AbstractSystem)
+
+Replicates the behaviour of `unknowns`, but ignores unknowns of subsystems.
+"""
+function unknowns_toplevel(sys::AbstractSystem)
+    if has_parent(sys) && (parent = get_parent(sys)) !== nothing
+        return unknowns_toplevel(parent)
+    end
+    return get_unknowns(sys)
+end
+
 """
 $(TYPEDSIGNATURES)
 
@@ -1353,6 +1365,18 @@ function dependent_parameters(sys::AbstractSystem)
     return map(eq -> eq.lhs, parameter_dependencies(sys))
 end
 
+"""
+    parameters_toplevel(sys::AbstractSystem)
+
+Replicates the behaviour of `parameters`, but ignores parameters of subsystems.
+"""
+function parameters_toplevel(sys::AbstractSystem)
+    if has_parent(sys) && (parent = get_parent(sys)) !== nothing
+        return parameters_toplevel(parent)
+    end
+    return get_ps(sys)
+end
+
 """
 $(TYPEDSIGNATURES)
 Get the parameter dependencies of the system `sys` and its subsystems.
@@ -1563,6 +1587,22 @@ function equations(sys::AbstractSystem)
     end
 end
 
+"""
+    equations_toplevel(sys::AbstractSystem)
+
+Replicates the behaviour of `equations`, but ignores equations of subsystems.
+
+Notes:
+- Cannot be applied to non-complete systems.
+"""
+function equations_toplevel(sys::AbstractSystem)
+    iscomplete(sys) && error("Cannot apply `equations_toplevel` to complete systems.")
+    if has_parent(sys) && (parent = get_parent(sys)) !== nothing
+        return equations_toplevel(parent)
+    end
+    return get_eqs(sys)
+end
+
 """
 $(TYPEDSIGNATURES)
 
diff --git a/src/systems/callbacks.jl b/src/systems/callbacks.jl
index 5cefe65651..937264d083 100644
--- a/src/systems/callbacks.jl
+++ b/src/systems/callbacks.jl
@@ -94,17 +94,17 @@ A [`ContinuousCallback`](@ref SciMLBase.ContinuousCallback) specified symbolical
 as well as the positive-edge `affect` and negative-edge `affect_neg` that apply when *any* of `eq` are satisfied.
 By default `affect_neg = affect`; to only get rising edges specify `affect_neg = nothing`.
 
-Assume without loss of generality that the equation is of the form `c(u,p,t) ~ 0`; we denote the integrator state as `i.u`. 
+Assume without loss of generality that the equation is of the form `c(u,p,t) ~ 0`; we denote the integrator state as `i.u`.
 For compactness, we define `prev_sign = sign(c(u[t-1], p[t-1], t-1))` and `cur_sign = sign(c(u[t], p[t], t))`.
-A condition edge will be detected and the callback will be invoked iff `prev_sign * cur_sign <= 0`. 
+A condition edge will be detected and the callback will be invoked iff `prev_sign * cur_sign <= 0`.
 The positive edge `affect` will be triggered iff an edge is detected and if `prev_sign < 0`; similarly, `affect_neg` will be
-triggered iff an edge is detected and `prev_sign > 0`. 
+triggered iff an edge is detected and `prev_sign > 0`.
 
-Inter-sample condition activation is not guaranteed; for example if we use the dirac delta function as `c` to insert a 
+Inter-sample condition activation is not guaranteed; for example if we use the dirac delta function as `c` to insert a
 sharp discontinuity between integrator steps (which in this example would not normally be identified by adaptivity) then the condition is not
 guaranteed to be triggered.
 
-Once detected the integrator will "wind back" through a root-finding process to identify the point when the condition became active; the method used 
+Once detected the integrator will "wind back" through a root-finding process to identify the point when the condition became active; the method used
 is specified by `rootfind` from [`SciMLBase.RootfindOpt`](@ref). If we denote the time when the condition becomes active as `tc`,
 the value in the integrator after windback will be:
 * `u[tc-epsilon], p[tc-epsilon], tc` if `LeftRootFind` is used,
@@ -116,7 +116,7 @@ it passed through 0.
 
 Multiple callbacks in the same system with different `rootfind` operations will be grouped
 by their `rootfind` value into separate VectorContinuousCallbacks in the enumeration order of `SciMLBase.RootfindOpt`. This may cause some callbacks to not fire if several become
-active at the same instant. See the `SciMLBase` documentation for more information on the semantic rules. 
+active at the same instant. See the `SciMLBase` documentation for more information on the semantic rules.
 
 Affects (i.e. `affect` and `affect_neg`) can be specified as either:
 * A list of equations that should be applied when the callback is triggered (e.g. `x ~ 3, y ~ 7`) which must be of the form `unknown ~ observed value` where each `unknown` appears only once. Equations will be applied in the order that they appear in the vector; parameters and state updates will become immediately visible to following equations.
@@ -328,13 +328,13 @@ The `SymbolicContinuousCallback`s in the returned vector are structs with two fi
 `eqs => affect`.
 """
 function continuous_events(sys::AbstractSystem)
-    obs = get_continuous_events(sys)
-    filter(!isempty, obs)
+    cbs = get_continuous_events(sys)
+    filter(!isempty, cbs)
 
     systems = get_systems(sys)
-    cbs = [obs;
+    cbs = [cbs;
            reduce(vcat,
-               (map(o -> namespace_callback(o, s), continuous_events(s))
+               (map(cb -> namespace_callback(cb, s), continuous_events(s))
                for s in systems),
                init = SymbolicContinuousCallback[])]
     filter(!isempty, cbs)
@@ -356,6 +356,21 @@ function vars!(vars, cb::SymbolicContinuousCallback; op = Differential)
     return vars
 end
 
+"""
+    continuous_events_toplevel(sys::AbstractSystem)
+
+Replicates the behaviour of `continuous_events`, but ignores events of subsystems.
+
+Notes:
+- Cannot be applied to non-complete systems.
+"""
+function continuous_events_toplevel(sys::AbstractSystem)
+    if has_parent(sys) && (parent = get_parent(sys)) !== nothing
+        return continuous_events_toplevel(parent)
+    end
+    return get_continuous_events(sys)
+end
+
 #################################### discrete events #####################################
 
 struct SymbolicDiscreteCallback
@@ -483,11 +498,11 @@ The `SymbolicDiscreteCallback`s in the returned vector are structs with two fiel
 `condition => affect`.
 """
 function discrete_events(sys::AbstractSystem)
-    obs = get_discrete_events(sys)
+    cbs = get_discrete_events(sys)
     systems = get_systems(sys)
-    cbs = [obs;
+    cbs = [cbs;
            reduce(vcat,
-               (map(o -> namespace_callback(o, s), discrete_events(s)) for s in systems),
+               (map(cb -> namespace_callback(cb, s), discrete_events(s)) for s in systems),
                init = SymbolicDiscreteCallback[])]
     cbs
 end
@@ -514,6 +529,21 @@ function vars!(vars, cb::SymbolicDiscreteCallback; op = Differential)
     return vars
 end
 
+"""
+    discrete_events_toplevel(sys::AbstractSystem)
+
+Replicates the behaviour of `discrete_events`, but ignores events of subsystems.
+
+Notes:
+- Cannot be applied to non-complete systems.
+"""
+function discrete_events_toplevel(sys::AbstractSystem)
+    if has_parent(sys) && (parent = get_parent(sys)) !== nothing
+        return discrete_events_toplevel(parent)
+    end
+    return get_discrete_events(sys)
+end
+
 ################################# compilation functions ####################################
 
 # handles ensuring that affect! functions work with integrator arguments
@@ -688,7 +718,7 @@ function generate_rootfinding_callback(sys::AbstractTimeDependentSystem,
     generate_rootfinding_callback(cbs, sys, dvs, ps; kwargs...)
 end
 """
-Generate a single rootfinding callback; this happens if there is only one equation in `cbs` passed to 
+Generate a single rootfinding callback; this happens if there is only one equation in `cbs` passed to
 generate_rootfinding_callback and thus we can produce a ContinuousCallback instead of a VectorContinuousCallback.
 """
 function generate_single_rootfinding_callback(
diff --git a/test/accessor_functions.jl b/test/accessor_functions.jl
new file mode 100644
index 0000000000..7ce477155b
--- /dev/null
+++ b/test/accessor_functions.jl
@@ -0,0 +1,171 @@
+### Preparations ###
+
+# Fetch packages.
+using ModelingToolkit, Test
+using ModelingToolkit: t_nounits as t, D_nounits as D
+import ModelingToolkit: get_ps, get_unknowns, get_observed, get_eqs, get_continuous_events,
+                        get_discrete_events, namespace_equations
+import ModelingToolkit: parameters_toplevel, unknowns_toplevel, equations_toplevel,
+                        continuous_events_toplevel, discrete_events_toplevel
+
+# Creates helper functions.
+function all_sets_equal(args...)
+    for arg in args[2:end]
+        issetequal(args[1], arg) || return false
+    end
+    return true
+end
+function sym_issubset(set1, set2)
+    for sym1 in set1
+        any(isequal(sym1, sym2) for sym2 in set2) || return false
+    end
+    return true
+end
+
+### Basic Tests ###
+
+# Checks `toplevel = false` argument for various accessors (currently only for `ODESystem`s).
+# Compares to `` version, and `get_` functions.
+# Checks  accessors for parameters, unknowns, equations, observables, and events.
+# Some tests looks funny (caused by the formatter).
+let
+    # Prepares model components.
+    @parameters p_top p_mid1 p_mid2 p_bot d
+    @variables X_top(t) X_mid1(t) X_mid2(t) X_bot(t) Y(t) O(t)
+
+    # Creates the systems (individual and hierarchical).
+    eqs_top = [
+        D(X_top) ~ p_top - d * X_top,
+        D(Y) ~ log(X_top) - Y^2 + 3.0,
+        O ~ (p_top + d) * X_top + Y
+    ]
+    eqs_mid1 = [
+        D(X_mid1) ~ p_mid1 - d * X_mid1^2,
+        D(Y) ~ D(X_mid1) - Y^3,
+        O ~ (p_mid1 + d) * X_mid1 + Y
+    ]
+    eqs_mid2 = [
+        D(X_mid2) ~ p_mid2 - d * X_mid2,
+        X_mid2^3 ~ log(X_mid2 + Y) - Y^2 + 3.0,
+        O ~ (p_mid2 + d) * X_mid2 + Y
+    ]
+    eqs_bot = [
+        D(X_bot) ~ p_bot - d * X_bot,
+        D(Y) ~ -Y^3,
+        O ~ (p_bot + d) * X_bot + Y
+    ]
+    cevs = [[t ~ 1.0] => [Y ~ Y + 2.0]]
+    devs = [(t == 2.0) => [Y ~ Y + 2.0]]
+    @named sys_bot = ODESystem(
+        eqs_bot, t; systems = [], continuous_events = cevs, discrete_events = devs)
+    @named sys_mid2 = ODESystem(
+        eqs_mid2, t; systems = [], continuous_events = cevs, discrete_events = devs)
+    @named sys_mid1 = ODESystem(
+        eqs_mid1, t; systems = [sys_bot], continuous_events = cevs, discrete_events = devs)
+    @named sys_top = ODESystem(eqs_top, t; systems = [sys_mid1, sys_mid2],
+        continuous_events = cevs, discrete_events = devs)
+    sys_bot_comp = complete(sys_bot)
+    sys_mid2_comp = complete(sys_mid2)
+    sys_mid1_comp = complete(sys_mid1)
+    sys_top_comp = complete(sys_top)
+    sys_bot_ss = structural_simplify(sys_bot)
+    sys_mid2_ss = structural_simplify(sys_mid2)
+    sys_mid1_ss = structural_simplify(sys_mid1)
+    sys_top_ss = structural_simplify(sys_top)
+
+    # Checks `parameters1.
+    @test all_sets_equal(parameters.([sys_bot, sys_bot_comp, sys_bot_ss])..., [d, p_bot])
+    @test all_sets_equal(parameters.([sys_mid1, sys_mid1_comp, sys_mid1_ss])...,
+        [d, p_mid1, sys_bot.d, sys_bot.p_bot])
+    @test all_sets_equal(
+        parameters.([sys_mid2, sys_mid2_comp, sys_mid2_ss])..., [d, p_mid2])
+    @test all_sets_equal(parameters.([sys_top, sys_top_comp, sys_top_ss])...,
+        [d, p_top, sys_mid1.d, sys_mid1.p_mid1, sys_mid1.sys_bot.d,
+            sys_mid1.sys_bot.p_bot, sys_mid2.d, sys_mid2.p_mid2])
+
+    # Checks `parameters_toplevel`. Compares to known parameters and also checks that
+    # these are subset of what `get_ps` returns.
+    @test all_sets_equal(
+        parameters_toplevel.([sys_bot, sys_bot_comp, sys_bot_ss])..., [d, p_bot])
+    @test all_sets_equal(
+        parameters_toplevel.([sys_mid1, sys_mid1_comp, sys_mid1_ss])...,
+        [d, p_mid1])
+    @test all_sets_equal(
+        parameters_toplevel.([sys_mid2, sys_mid2_comp, sys_mid2_ss])...,
+        [d, p_mid2])
+    @test all_sets_equal(
+        parameters_toplevel.([sys_top, sys_top_comp, sys_top_ss])..., [d, p_top])
+    @test all(sym_issubset(parameters_toplevel(sys), get_ps(sys))
+    for sys in [sys_bot, sys_mid2, sys_mid1, sys_top])
+
+    # Checks `unknowns`. O(t) is eliminated by `structural_simplify` and
+    # must be considered separately.
+    @test all_sets_equal(unknowns.([sys_bot, sys_bot_comp])..., [O, Y, X_bot])
+    @test all_sets_equal(unknowns.([sys_bot_ss])..., [Y, X_bot])
+    @test all_sets_equal(unknowns.([sys_mid1, sys_mid1_comp])...,
+        [O, Y, X_mid1, sys_bot.Y, sys_bot.O, sys_bot.X_bot])
+    @test all_sets_equal(unknowns.([sys_mid1_ss])..., [Y, X_mid1, sys_bot.Y, sys_bot.X_bot])
+    @test all_sets_equal(unknowns.([sys_mid2, sys_mid2_comp])..., [O, Y, X_mid2])
+    @test all_sets_equal(unknowns.([sys_mid2_ss])..., [Y, X_mid2])
+    @test all_sets_equal(unknowns.([sys_top, sys_top_comp])...,
+        [O, Y, X_top, sys_mid1.O, sys_mid1.Y, sys_mid1.X_mid1,
+            sys_mid1.sys_bot.O, sys_mid1.sys_bot.Y, sys_mid1.sys_bot.X_bot,
+            sys_mid2.O, sys_mid2.Y, sys_mid2.X_mid2])
+    @test all_sets_equal(unknowns.([sys_top_ss])...,
+        [Y, X_top, sys_mid1.Y, sys_mid1.X_mid1, sys_mid1.sys_bot.Y,
+            sys_mid1.sys_bot.X_bot, sys_mid2.Y, sys_mid2.X_mid2])
+
+    # Checks `unknowns_toplevel`. Note that O is not eliminated here (as we go back
+    # to original parent system). Also checks that outputs are subsets of what `get_unknowns` returns.
+    @test all_sets_equal(
+        unknowns_toplevel.([sys_bot, sys_bot_comp, sys_bot_ss])..., [O, Y, X_bot])
+    @test all_sets_equal(
+        unknowns_toplevel.([sys_mid1, sys_mid1_comp])..., [O, Y, X_mid1])
+    @test all_sets_equal(
+        unknowns_toplevel.([sys_mid2, sys_mid2_comp])..., [O, Y, X_mid2])
+    @test all_sets_equal(
+        unknowns_toplevel.([sys_top, sys_top_comp])..., [O, Y, X_top])
+    @test all(sym_issubset(unknowns_toplevel(sys), get_unknowns(sys))
+    for sys in [sys_bot, sys_mid1, sys_mid2, sys_top])
+
+    # Checks `equations`. Do not check ss equations as these might potentially
+    # be structurally simplified to new equations.
+    @test all_sets_equal(equations.([sys_bot, sys_bot_comp])..., eqs_bot)
+    @test all_sets_equal(
+        equations.([sys_mid1, sys_mid1_comp])..., [eqs_mid1; namespace_equations(sys_bot)])
+    @test all_sets_equal(equations.([sys_mid2, sys_mid2_comp])..., eqs_mid2)
+    @test all_sets_equal(equations.([sys_top, sys_top_comp])...,
+        [eqs_top; namespace_equations(sys_mid1); namespace_equations(sys_mid2)])
+
+    # Checks `equations_toplevel`. Do not check ss equations directly as these
+    # might potentially be structurally simplified to new equations. Do not check
+    @test all_sets_equal(equations_toplevel.([sys_bot])..., eqs_bot)
+    @test all_sets_equal(
+        equations_toplevel.([sys_mid1])..., eqs_mid1)
+    @test all_sets_equal(
+        equations_toplevel.([sys_mid2])..., eqs_mid2)
+    @test all_sets_equal(equations_toplevel.([sys_top])..., eqs_top)
+    @test all(sym_issubset(equations_toplevel(sys), get_eqs(sys))
+    for sys in [sys_bot, sys_mid2, sys_mid1, sys_top])
+
+    # Checks `continuous_events_toplevel` and  `discrete_events_toplevel` (straightforward
+    # as I stored the same singe event in all systems). Don't check for non-toplevel cases as
+    # technically not needed for these tests and name spacing the events is a mess.
+    mtk_cev = ModelingToolkit.SymbolicContinuousCallback.(cevs)[1]
+    mtk_dev = ModelingToolkit.SymbolicDiscreteCallback.(devs)[1]
+    @test all_sets_equal(
+        continuous_events_toplevel.(
+            [sys_bot, sys_bot_comp, sys_bot_ss, sys_mid1, sys_mid1_comp, sys_mid1_ss,
+            sys_mid2, sys_mid2_comp, sys_mid2_ss, sys_top, sys_top_comp, sys_top_ss])...,
+        [mtk_cev])
+    @test all_sets_equal(
+        discrete_events_toplevel.(
+            [sys_bot, sys_bot_comp, sys_bot_ss, sys_mid1, sys_mid1_comp, sys_mid1_ss,
+            sys_mid2, sys_mid2_comp, sys_mid2_ss, sys_top, sys_top_comp, sys_top_ss])...,
+        [mtk_dev])
+    @test all(sym_issubset(
+                  continuous_events_toplevel(sys), get_continuous_events(sys))
+    for sys in [sys_bot, sys_mid2, sys_mid1, sys_top])
+    @test all(sym_issubset(discrete_events_toplevel(sys), get_discrete_events(sys))
+    for sys in [sys_bot, sys_mid2, sys_mid1, sys_top])
+end
diff --git a/test/runtests.jl b/test/runtests.jl
index e0ef4b8640..4537c27255 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -61,6 +61,7 @@ end
             @safetestset "Constants Test" include("constants.jl")
             @safetestset "Parameter Dependency Test" include("parameter_dependencies.jl")
             @safetestset "Equation Type Accessors Test" include("equation_type_accessors.jl")
+            @safetestset "System Accessor Functions Test" include("accessor_functions.jl")
             @safetestset "Equations with complex values" include("complex.jl")
         end
     end