diff --git a/.github/workflows/ci.yml b/.github/workflows/ci.yml
index 0fb0a7b3..246886c4 100644
--- a/.github/workflows/ci.yml
+++ b/.github/workflows/ci.yml
@@ -21,7 +21,7 @@ jobs:
       fail-fast: false
       matrix:
         julia-version:
-          - "1.0"
+          - "1.6"
           - "1"
           - "nightly"
         os:
diff --git a/NEWS.md b/NEWS.md
index 80e5b448..fe4674bb 100644
--- a/NEWS.md
+++ b/NEWS.md
@@ -1,5 +1,12 @@
 # History of Measurements.jl
 
+## v2.12.0 (????)
+
+### New features
+
+* Support for [symbolic variables](https://symbolics.juliasymbolics.org/stable/manual/variables/),
+   as provided by `Symbolics.jl` package.
+
 ## v2.11.0 (2023-11-03)
 
 ### New features
diff --git a/Project.toml b/Project.toml
index e3ff5c49..91e2f298 100644
--- a/Project.toml
+++ b/Project.toml
@@ -15,6 +15,7 @@ BaseType = "7fbed51b-1ef5-4d67-9085-a4a9b26f478c"
 Juno = "e5e0dc1b-0480-54bc-9374-aad01c23163d"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
+Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [extensions]
@@ -22,6 +23,7 @@ MeasurementsBaseTypeExt = "BaseType"
 MeasurementsJunoExt = "Juno"
 MeasurementsRecipesBaseExt = "RecipesBase"
 MeasurementsSpecialFunctionsExt = "SpecialFunctions"
+MeasurementsSymbolicsExt = "Symbolics"
 MeasurementsUnitfulExt = "Unitful"
 
 [compat]
@@ -30,7 +32,7 @@ Calculus = "0.4.1, 0.5"
 LinearAlgebra = "<0.0.1, 1"
 RecipesBase = "0.6.0, 0.7, 0.8, 1.0"
 Requires = "0.5.0, 1"
-julia = "1"
+julia = "1.6"
 
 [extras]
 Aqua = "4c88cf16-eb10-579e-8560-4a9242c79595"
@@ -40,8 +42,9 @@ QuadGK = "1fd47b50-473d-5c70-9696-f719f8f3bcdc"
 RecipesBase = "3cdcf5f2-1ef4-517c-9805-6587b60abb01"
 SpecialFunctions = "276daf66-3868-5448-9aa4-cd146d93841b"
 Statistics = "10745b16-79ce-11e8-11f9-7d13ad32a3b2"
+Symbolics = "0c5d862f-8b57-4792-8d23-62f2024744c7"
 Test = "8dfed614-e22c-5e08-85e1-65c5234f0b40"
 Unitful = "1986cc42-f94f-5a68-af5c-568840ba703d"
 
 [targets]
-test = ["Aqua", "BaseType", "QuadGK", "RecipesBase", "SpecialFunctions", "Statistics", "Test", "Unitful"]
+test = ["Aqua", "BaseType", "QuadGK", "RecipesBase", "SpecialFunctions", "Statistics", "Symbolics", "Test", "Unitful"]
diff --git a/README.md b/README.md
index 92b05fbc..40a7c00b 100644
--- a/README.md
+++ b/README.md
@@ -52,6 +52,9 @@ easy-to-use calculator.
   [arbitrary precision](https://docs.julialang.org/en/v1/manual/integers-and-floating-point-numbers/#Arbitrary-Precision-Arithmetic-1)
   (also called multiple precision) numbers with uncertainties.  This is useful
   for measurements with very low relative error
+* Support for numbers with uncertainties of
+  [symbolic variable](https://symbolics.juliasymbolics.org/stable/manual/variables/)
+  type, as provided by `Symbolics.jl` package.
 * Define arrays of measurements and perform calculations with them.  Some linear
   algebra functions work out-of-the-box
 * Propagate uncertainty for any function of real arguments (including functions
diff --git a/docs/src/appendix.md b/docs/src/appendix.md
index 2a078fc4..04d2b119 100644
--- a/docs/src/appendix.md
+++ b/docs/src/appendix.md
@@ -11,7 +11,7 @@ The `Measurement` Type
 `Measurement` is a
 [composite](https://docs.julialang.org/en/v1/manual/types/#Composite-Types-1)
 [parametric](https://docs.julialang.org/en/v1/manual/types/#Parametric-Types-1)
-type, whose parameter is the `AbstractFloat` subtype of the nominal value and
+type, whose parameter is the `Real` subtype of the nominal value and
 the uncertainty of the measurement. `Measurement` type itself is subtype of
 `AbstractFloat`, thus `Measurement` objects can be used in any function taking
 `AbstractFloat` arguments without redefining it, and calculation of uncertainty
@@ -20,7 +20,7 @@ will be exact.
 In detail, this is the definition of the type:
 
 ```julia
-struct Measurement{T<:AbstractFloat} <: AbstractFloat
+struct Measurement{T<:Real} <: AbstractFloat
     val::T
     err::T
     tag::UInt64
diff --git a/docs/src/examples.md b/docs/src/examples.md
index 3dd0633b..8ea41d2b 100644
--- a/docs/src/examples.md
+++ b/docs/src/examples.md
@@ -272,6 +272,33 @@ hypot(a, b)
 log(2a) ^ b
 ```
 
+Symbolic Variables
+------------------
+
+You can perform any mathematical operation supported in `Measurements.jl` using
+[symbolic variables](https://symbolics.juliasymbolics.org/stable/manual/variables/),
+as provided by `Symbolics.jl` package:
+
+``` julia
+julia> using Measurements, Symbolics
+
+julia> @variables x_val, x_err, y_val, y_err
+4-element Vector{Num}:
+ x_val x_err y_val y_err
+
+julia> x = x_val ± x_err
+x_val ± x_err
+
+julia> y = y_val ± y_err
+y_val ± y_err
+
+julia> x + y
+x_val + y_val ± sqrt(abs2(x_err) + abs2(y_err))
+
+julia> x - x
+0 ± 0.0
+```
+
 Operations with Arrays and Linear Algebra
 -----------------------------------------
 
diff --git a/docs/src/index.md b/docs/src/index.md
index 9d881458..cb09a2b4 100644
--- a/docs/src/index.md
+++ b/docs/src/index.md
@@ -48,6 +48,9 @@ The main features of the package are:
   precision](https://docs.julialang.org/en/v1/manual/integers-and-floating-point-numbers/#Arbitrary-Precision-Arithmetic-1)
   (also called multiple precision) numbers with uncertainties. This is useful
   for measurements with very low relative error
+- Support for numbers with uncertainties of
+  [symbolic variable](https://symbolics.juliasymbolics.org/stable/manual/variables/)
+  type, as provided by `Symbolics.jl` package.
 - Define arrays of measurements and perform calculations with them.  Some linear
   algebra functions work out-of-the-box
 - Propagate uncertainty for any function of real arguments (including functions
diff --git a/docs/src/usage.md b/docs/src/usage.md
index 92265f08..3d541c79 100644
--- a/docs/src/usage.md
+++ b/docs/src/usage.md
@@ -20,11 +20,12 @@ where
 - `err` is its uncertainty, assumed to be a [standard
   deviation](https://en.wikipedia.org/wiki/Standard_deviation).
 
-They are both subtype of `AbstractFloat`. Some keyboard layouts provide an easy
+They are both subtype of `Real`. Some keyboard layouts provide an easy
 way to type the `±` sign, if your does not, remember you can insert it in Julia
 REPL with `\pm` followed by `TAB` key. You can provide `val` and `err` of any
-subtype of `Real` that can be converted to `AbstractFloat`. Thus,
-`measurement(42, 33//12)` and `pi ± 0.1` are valid.
+subtype of `Real` that can be converted to `AbstractFloat`, or are
+[symbolic variables](https://symbolics.juliasymbolics.org/stable/manual/variables/).
+Thus, `measurement(42, 33//12)` and `pi ± 0.1` are valid.
 
 `measurement(value)` creates a `Measurement` object with zero uncertainty, like
 mathematical constants. See below for further examples.
diff --git a/ext/MeasurementsSpecialFunctionsExt.jl b/ext/MeasurementsSpecialFunctionsExt.jl
index 8879f019..dbfb7c26 100644
--- a/ext/MeasurementsSpecialFunctionsExt.jl
+++ b/ext/MeasurementsSpecialFunctionsExt.jl
@@ -30,40 +30,40 @@ else
 end
 # Error function: erf, erfinv, erfc, erfcinv, erfcx, erfi, dawson
 
-function SpecialFunctions.erf(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erf(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(erf(aval), 2*exp(-abs2(aval))/sqrt(T(pi)), a)
 end
 
-function SpecialFunctions.erfinv(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erfinv(a::Measurement{T}) where {T<:Real}
     res = erfinv(a.val)
     # For the derivative, see http://mathworld.wolfram.com/InverseErf.html
     return result(res, sqrt(T(pi)) * exp(abs2(res)) / 2, a)
 end
 
-function SpecialFunctions.erfc(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erfc(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(erfc(aval), -2*exp(-abs2(aval))/sqrt(T(pi)), a)
 end
 
-function SpecialFunctions.erfcinv(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erfcinv(a::Measurement{T}) where {T<:Real}
     res = erfcinv(a.val)
     # For the derivative, see http://mathworld.wolfram.com/InverseErfc.html
     return result(res, -sqrt(T(pi)) * exp(abs2(res)) / 2, a)
 end
 
-function SpecialFunctions.erfcx(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erfcx(a::Measurement{T}) where {T<:Real}
     aval = a.val
     res = erfcx(aval)
     return result(res, 2 * (aval * res - inv(sqrt(T(pi)))), a)
 end
 
-function SpecialFunctions.erfi(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.erfi(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(erfi(aval), 2*exp(abs2(aval))/sqrt(T(pi)), a)
 end
 
-function SpecialFunctions.dawson(a::Measurement{T}) where {T<:AbstractFloat}
+function SpecialFunctions.dawson(a::Measurement{T}) where {T<:Real}
     aval = a.val
     res = dawson(aval)
     return result(res, one(T) - 2 * aval * res, a)
diff --git a/ext/MeasurementsSymbolicsExt.jl b/ext/MeasurementsSymbolicsExt.jl
new file mode 100644
index 00000000..2911840f
--- /dev/null
+++ b/ext/MeasurementsSymbolicsExt.jl
@@ -0,0 +1,30 @@
+### MeasurementsSymbolicsExt.jl
+#
+# Copyright (C) 2023 Roman Lebedev.
+#
+# Maintainer: Mosè Giordano <mose AT gnu DOT org>
+# Keywords: uncertainty, error propagation, physics
+#
+# This file is a part of Measurements.jl.
+#
+# License is MIT "Expat".
+#
+### Commentary:
+#
+# This file contains integration with Symbolics.jl.
+#
+### Code:
+
+### Symbolics
+module MeasurementsSymbolicsExt
+
+import Measurements
+import Symbolics
+
+Measurements._is_symbolic(::Symbolics.Num) = true
+
+Measurements.measurement(val::Symbolics.Num) = Measurements._measurement(val)
+
+Measurements.measurement(val::Symbolics.Num, err::Symbolics.Num) = Measurements._measurement(val, err)
+
+end
diff --git a/src/Measurements.jl b/src/Measurements.jl
index c75ed50d..2f40ce51 100644
--- a/src/Measurements.jl
+++ b/src/Measurements.jl
@@ -27,6 +27,9 @@ using Calculus
 # Functions provided by this package and exposed to users
 export Measurement, measurement, ±
 
+# Query whether the value is of symbolic type.
+_is_symbolic(::Real) = false
+
 # Define the "Derivatives" type, used inside "Measurement" type.  This should be
 # a lightweight and immutable dictionary.
 include("derivatives-type.jl")
@@ -46,7 +49,7 @@ include("derivatives-type.jl")
 #     measurement and propagate the uncertainty in the case of functions with
 #     more than one argument (in order to deal with correlation between
 #     arguments).
-struct Measurement{T<:AbstractFloat} <: AbstractFloat
+struct Measurement{T<:Real} <: AbstractFloat
     val::T
     err::T
     tag::UInt64
@@ -72,16 +75,19 @@ function Measurement{T}(::S) where {T, S}
 end
 
 # Functions to quickly create an empty Derivatives object.
-@generated empty_der1(x::Measurement{T}) where {T<:AbstractFloat} = Derivatives{T}()
-@generated empty_der2(x::T) where {T<:AbstractFloat} = Derivatives{x}()
+@generated empty_der1(x::Measurement{T}) where {T<:Real} = Derivatives{T}()
+@generated empty_der2(x::T) where {T<:Real} = Derivatives{x}()
 
 # Start from 1, 0 is reserved to derived quantities
 const tag_counter = Threads.Atomic{UInt64}(1)
 
 measurement(x::Measurement) = x
-measurement(val::T) where {T<:AbstractFloat} = Measurement(val, zero(T), UInt64(0), empty_der2(val))
+
+_measurement(val::T) where {T<:Real} =  Measurement(val, zero(T), UInt64(0), empty_der2(val))
+measurement(val::AbstractFloat) = _measurement(val)
 measurement(val::Real) = measurement(float(val))
-function measurement(val::T, err::T) where {T<:AbstractFloat}
+
+function _measurement(val::T, err::T) where {T<:Real}
     newder = empty_der2(val)
     if iszero(err)
         Measurement{T}(val, err, UInt64(0), newder)
@@ -90,8 +96,12 @@ function measurement(val::T, err::T) where {T<:AbstractFloat}
         return Measurement{T}(val, err, tag, Derivatives(newder, (val, err, tag)=>one(T)))
     end
 end
-measurement(val::Real, err::Real) = measurement(promote(float(val), float(err))...)
+measurement(val::T, err::T) where {T<:AbstractFloat} = _measurement(val, err)
+measurement(val::T, err::T) where {T<:Real} = measurement(float(val), float(err))
+measurement(val::V, err::E) where {V<:Real, E<:Real} = measurement(promote(val, err)...)
+
 measurement(::Missing, ::Union{Real,Missing} = missing) = missing
+
 const ± = measurement
 
 """
diff --git a/src/conversions.jl b/src/conversions.jl
index 753ebfbe..b1b6ac96 100644
--- a/src/conversions.jl
+++ b/src/conversions.jl
@@ -16,15 +16,15 @@
 #
 ### Code:
 
-Base.convert(::Type{Measurement{T}}, a::Irrational) where {T<:AbstractFloat} =
+Base.convert(::Type{Measurement{T}}, a::Irrational) where {T<:Real} =
     measurement(T(a))::Measurement{T}
-Base.convert(::Type{Measurement{T}}, a::Rational{<:Integer}) where {T<:AbstractFloat} =
+Base.convert(::Type{Measurement{T}}, a::Rational{<:Integer}) where {T<:Real} =
     measurement(T(a))::Measurement{T}
-Base.convert(::Type{Measurement{T}}, a::Real) where {T<:AbstractFloat} =
+Base.convert(::Type{Measurement{T}}, a::Real) where {T<:Real} =
     measurement(T(a))::Measurement{T}
-Base.convert(::Type{Measurement{T}}, a::Base.TwicePrecision) where {T<:AbstractFloat} =
+Base.convert(::Type{Measurement{T}}, a::Base.TwicePrecision) where {T<:Real} =
     measurement(T(a))::Measurement{T}
-Base.convert(::Type{Measurement{T}}, a::AbstractChar) where {T<:AbstractFloat} =
+Base.convert(::Type{Measurement{T}}, a::AbstractChar) where {T<:Real} =
     measurement(T(a))::Measurement{T}
 
 function Base.convert(::Type{Measurement{T}}, a::Complex) where {T}
@@ -35,9 +35,9 @@ function Base.convert(::Type{Measurement{T}}, a::Complex) where {T}
     end
 end
 
-Base.convert(::Type{Measurement{T}}, a::Measurement{T}) where {T<:AbstractFloat} = a
+Base.convert(::Type{Measurement{T}}, a::Measurement{T}) where {T<:Real} = a
 function Base.convert(::Type{Measurement{T}},
-                      a::Measurement{<:AbstractFloat}) where {T<:AbstractFloat}
+                      a::Measurement{<:Real}) where {T<:Real}
     newder = empty_der2(zero(T))
     for tag in keys(a.der)
         newder = Derivatives(newder, (T(tag[1]), T(tag[2]), tag[3])=>T(a.der[tag]))
@@ -55,10 +55,10 @@ function Base.convert(::Type{Int}, a::Measurement)
     return convert(Int, a.val)::Int
 end
 
-Base.promote_rule(::Type{Measurement{T}}, ::Type{S}) where {T<:AbstractFloat, S<:Real} =
+Base.promote_rule(::Type{Measurement{T}}, ::Type{S}) where {T<:Real, S<:Real} =
     Measurement{promote_type(T, S)}
 Base.promote_rule(::Type{Measurement{T}},
-                  ::Type{Measurement{S}}) where {T<:AbstractFloat, S<:AbstractFloat} =
+                  ::Type{Measurement{S}}) where {T<:Real, S<:Real} =
     Measurement{promote_type(T, S)}
 
 # adaptation of JuliaLang/julia#30952
diff --git a/src/math.jl b/src/math.jl
index 3115e7f1..6e97b341 100644
--- a/src/math.jl
+++ b/src/math.jl
@@ -38,7 +38,7 @@ export @uncertain
 #   σ_G = |σ_a·∂G/∂a|
 # The list of derivatives with respect to each measurement is updated with
 #   ∂G/∂a · previous_derivatives
-@inline function result(val::T, der::Real, a::Measurement{<:AbstractFloat}) where {T<:Real}
+@inline function result(val::T, der::Real, a::Measurement{<:Real}) where {T<:Real}
     newder = empty_der1(a)
     @inbounds for tag in keys(a.der)
         if ! iszero(tag[2]) # Skip values with 0 uncertainty
@@ -56,9 +56,9 @@ end
 # Get the common type parameter of a collection of Measurement objects.  The first two
 # methods are for the trivial cases of homogeneous tuples and arrays, the last, inefficient,
 # method is for inhomogeneous collections (probably the least common case).
-gettype(::Tuple{Measurement{T}, Vararg{Measurement{T}}}) where {T<:AbstractFloat} = T
-gettype(::AbstractArray{Measurement{T}}) where {T<:AbstractFloat} = T
-_eltype(::Measurement{T}) where {T<:AbstractFloat} = T
+gettype(::Tuple{Measurement{T}, Vararg{Measurement{T}}}) where {T<:Real} = T
+gettype(::AbstractArray{Measurement{T}}) where {T<:Real} = T
+_eltype(::Measurement{T}) where {T<:Real} = T
 gettype(collection) = promote_type(_eltype.(collection)...)
 
 # This function is similar to the previous one, but applies to mathematical
@@ -249,7 +249,7 @@ function Base.:/(a::Measurement, b::Measurement)
     return result(x / y, (oneovery, -x * abs2(oneovery)), (a, b))
 end
 Base.:/(a::Real, b::Measurement) = result(a/b.val, -a/abs2(b.val), b)
-Base.:/(a::Measurement{T}, b::Real) where {T<:AbstractFloat} = result(a.val/b, 1/T(b), a)
+Base.:/(a::Measurement{T}, b::Real) where {T<:Real} = result(a.val/b, 1/T(b), a)
 
 # 0.0 as partial derivative for both arguments of "div", "fld", "cld" should be
 # correct for most cases.  This has been tested against "@uncertain" macro.
@@ -289,7 +289,7 @@ function Base.:^(a::Measurement, b::Integer)
     return result(x ^ b, b * x ^ (b - 1), a)
 end
 
-function Base.:^(a::Measurement{T},  b::Rational) where {T<:AbstractFloat}
+function Base.:^(a::Measurement{T},  b::Rational) where {T<:Real}
     x = a.val
     return result(x ^ b, b * x ^ (b - one(T)), a)
 end
@@ -306,7 +306,7 @@ function Base.:^(a::Real,  b::Measurement)
     return result(res, res*log(a), b)
 end
 
-function Base.exp2(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.exp2(a::Measurement{T}) where {T<:Real}
     pow = exp2(a.val)
     return result(pow, pow*log(T(2)), a)
 end
@@ -411,47 +411,47 @@ end
 # Inverse trig functions: acos, acosd, acosh, asin, asind, asinh, atan, atand, atanh,
 #                         asec, acsc, acot, asech, acsch, acoth
 
-function Base.acos(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.acos(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(acos(aval), -inv(sqrt(one(T) - abs2(aval))), a)
 end
 
-function Base.acosd(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.acosd(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(acosd(aval), -rad2deg(inv(sqrt(one(T) - abs2(aval)))), a)
 end
 
-function Base.acosh(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.acosh(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(acosh(aval), inv(sqrt(abs2(aval) - one(T))), a)
 end
 
-function Base.asin(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.asin(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(asin(aval), inv(sqrt(one(T) - abs2(aval))), a)
 end
 
-function Base.asind(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.asind(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(asind(aval), rad2deg(inv(sqrt(one(T) - abs2(aval)))), a)
 end
 
-function Base.asinh(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.asinh(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(asinh(aval), inv(hypot(aval, one(T))), a)
 end
 
-function Base.atan(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.atan(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(atan(aval), inv(abs2(aval) + one(T)), a)
 end
 
-function Base.atand(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.atand(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(atand(aval), rad2deg(inv(abs2(aval) + one(T))), a)
 end
 
-function Base.atanh(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.atanh(a::Measurement{T}) where {T<:Real}
     aval = a.val
     return result(atanh(aval), inv(one(T) - abs2(aval)), a)
 end
@@ -572,7 +572,7 @@ function Base.expm1(a::Measurement)
     return result(expm1(aval), exp(aval), a)
 end
 
-function Base.exp10(a::Measurement{T}) where {T<:AbstractFloat}
+function Base.exp10(a::Measurement{T}) where {T<:Real}
     val = exp10(a.val)
     return result(val, log(T(10))*val, a)
 end
@@ -582,7 +582,7 @@ function Base.frexp(a::Measurement)
     return (result(x, inv(exp2(y)), a), y)
 end
 
-Base.ldexp(a::Measurement{T}, e::Integer) where {T<:AbstractFloat} =
+Base.ldexp(a::Measurement{T}, e::Integer) where {T<:Real} =
     result(ldexp(a.val, e), ldexp(one(T), e), a)
 
 # Logarithms
@@ -600,17 +600,17 @@ function Base.log(a::Measurement) # Special case
     return result(log(aval), inv(aval), a)
 end
 
-function Base.log2(a::Measurement{T}) where {T<:AbstractFloat} # Special case
+function Base.log2(a::Measurement{T}) where {T<:Real} # Special case
     x = a.val
     return result(log2(x), inv(log(T(2)) * x), a)
 end
 
-function Base.log10(a::Measurement{T}) where {T<:AbstractFloat} # Special case
+function Base.log10(a::Measurement{T}) where {T<:Real} # Special case
     aval = a.val
     return result(log10(aval), inv(log(T(10)) * aval), a)
 end
 
-function Base.log1p(a::Measurement{T}) where {T<:AbstractFloat} # Special case
+function Base.log1p(a::Measurement{T}) where {T<:Real} # Special case
     aval = a.val
     return result(log1p(aval), inv(aval + one(T)), a)
 end
@@ -719,18 +719,18 @@ Base.rem2pi(a::Measurement, r::RoundingMode) = result(rem2pi(a.val, r), 1, a)
 
 ### Machine precision
 
-Base.eps(::Type{Measurement{T}}) where {T<:AbstractFloat} = eps(T)
+Base.eps(::Type{Measurement{T}}) where {T<:Real} = eps(T)
 Base.eps(a::Measurement) = eps(a.val)
 
 Base.nextfloat(a::Measurement) = result(nextfloat(a.val), 1, a)
 Base.nextfloat(a::Measurement, n::Integer) = result(nextfloat(a.val, n), 1, a)
 
-Base.maxintfloat(::Type{Measurement{T}}) where {T<:AbstractFloat} = maxintfloat(T)
+Base.maxintfloat(::Type{Measurement{T}}) where {T<:Real} = maxintfloat(T)
 
-Base.floatmin(::Type{Measurement{T}}) where {T<:AbstractFloat} = floatmin(T) ± zero(T)
-Base.floatmax(::Type{Measurement{T}}) where {T<:AbstractFloat} = floatmax(T) ± zero(T)
+Base.floatmin(::Type{Measurement{T}}) where {T<:Real} = floatmin(T) ± zero(T)
+Base.floatmax(::Type{Measurement{T}}) where {T<:Real} = floatmax(T) ± zero(T)
 
-Base.typemax(::Type{Measurement{T}}) where {T<:AbstractFloat} = typemax(T)
+Base.typemax(::Type{Measurement{T}}) where {T<:Real} = typemax(T)
 
 ### Rounding
 
@@ -766,13 +766,13 @@ Base.trunc(::Type{Bool}, x::Measurement) = measurement(trunc(Bool, value(x)))
 
 # Widening
 
-Base.widen(::Type{Measurement{T}}) where {T<:AbstractFloat} = Measurement{widen(T)}
+Base.widen(::Type{Measurement{T}}) where {T<:Real} = Measurement{widen(T)}
 
 # To big float
 
 Base.big(::Type{Measurement}) = Measurement{BigFloat}
-Base.big(::Type{Measurement{T}}) where {T<:AbstractFloat} = Measurement{BigFloat}
-Base.big(x::Measurement{<:AbstractFloat}) = convert(Measurement{BigFloat}, x)
+Base.big(::Type{Measurement{T}}) where {T<:Real} = Measurement{BigFloat}
+Base.big(x::Measurement{<:Real}) = convert(Measurement{BigFloat}, x)
 Base.big(x::Complex{<:Measurement}) = convert(Complex{Measurement{BigFloat}}, x)
 
 # Sum and prod
diff --git a/src/parsing.jl b/src/parsing.jl
index 0488409b..6de9ea7c 100644
--- a/src/parsing.jl
+++ b/src/parsing.jl
@@ -83,7 +83,7 @@ julia> measurement("-1234e-1")
 """
 measurement(str::AbstractString) = parse(Measurement{Float64}, str)
 
-function Base.tryparse(::Type{Measurement{T}}, str::S) where {T<:AbstractFloat, S<:AbstractString}
+function Base.tryparse(::Type{Measurement{T}}, str::S) where {T<:Real, S<:AbstractString}
     m = match(rxp_error_with_parentheses, str)
     if m !== nothing # "123(45)e6"
         val_str::S, val_dec, err_str::S, err_dec_str, expn = m.captures
@@ -128,7 +128,7 @@ function Base.tryparse(::Type{Measurement{T}}, str::S) where {T<:AbstractFloat,
     return measurement(val, err)
 end
 
-function Base.parse(::Type{Measurement{T}}, str::S) where {T<:AbstractFloat, S<:AbstractString}
+function Base.parse(::Type{Measurement{T}}, str::S) where {T<:Real, S<:AbstractString}
     out = tryparse(Measurement{T}, str)
     out === nothing && throw(ArgumentError("cannot parse $(repr(str)) as Measurement{$T}"))
     return out
diff --git a/src/show.jl b/src/show.jl
index 033e7b46..241f747b 100644
--- a/src/show.jl
+++ b/src/show.jl
@@ -17,47 +17,46 @@
 
 import Printf
 
-function truncated_print(io::IO, m::Measurement, error_digits::Int;
+full_print(io::IO, val, err; atbeg = "", atend = "", pm = "±") =
+    print(io, atbeg, val, get(io, :compact, false) ? pm : " $pm ", err, atend)
+
+function truncated_print(io::IO, m::Measurement;
                          atbeg = "", atend = "", pm = "±")
-    val = if iszero(m.err) || !isfinite(m.err)
-        m.val
-    else
-        err_digits = -Base.hidigit(m.err, 10) + error_digits
-        digits = if isfinite(m.val)
-            max(-Base.hidigit(m.val, 10) + 2, err_digits)
+    error_digits::Int = get(io, :error_digits, 2)
+
+    val = m.val
+    err = m.err
+
+    if error_digits < 0
+        error("`error_digits` must be non-negative")
+    elseif error_digits > 0 && !_is_symbolic(val)
+        val = if iszero(err) || !isfinite(err)
+            val
         else
-            err_digits
+            err_digits = -Base.hidigit(err, 10) + error_digits
+            digits = if isfinite(val)
+                max(-Base.hidigit(val, 10) + 2, err_digits)
+            else
+                err_digits
+            end
+            round(val, digits = digits)
         end
-        round(m.val, digits = digits)
-    end
-    print(io, atbeg, val,
-          get(io, :compact, false) ? pm : " $pm ",
-          round(m.err, sigdigits=error_digits), atend)
+        err = round(err, sigdigits=error_digits)
+    end # if
+    full_print(io, val, err, atbeg = atbeg, atend = atend, pm = pm)
 end
 
-full_print(io::IO, measure::Measurement) =
-    print(io, measure.val, get(io, :compact, false) ? "±" : " ± ", measure.err)
-
 function Base.show(io::IO, m::Measurement)
-    error_digits = get(io, :error_digits, 2)
-    if error_digits > 0
-            truncated_print(io, m, error_digits)
-    elseif error_digits == 0
-        full_print(io, m)
-    else
-        error("`error_digits` must be non-negative")
-    end # if
+    truncated_print(io, m)
 end
 
 function Base.show(io::IO, ::MIME"text/latex", measure::Measurement)
-    error_digits = get(io, :error_digits, 2)
-    truncated_print(io, measure, error_digits, atbeg = "\$", atend = "\$", pm = "\\pm")
+    truncated_print(io, measure, atbeg = "\$", atend = "\$", pm = "\\pm")
 end
 
 for mime in (MIME"text/x-tex", MIME"text/x-latex")
     function Base.show(io::IO, ::mime, measure::Measurement)
-        error_digits = get(io, :error_digits, 2)
-        truncated_print(io, measure, error_digits, pm = "\\pm")
+        truncated_print(io, measure, pm = "\\pm")
     end
 end
 
@@ -71,7 +70,7 @@ for mime in (MIME"text/plain", MIME"text/x-tex", MIME"text/x-latex")
         print(io, "(")
         show(io, mtype, r)
         print(io, ")")
-        if signbit(i) && !isnan(i)
+        if !_is_symbolic(i.val) && signbit(i) && !isnan(i)
             i = -i
             print(io, compact ? "-" : " - ")
         else
@@ -83,6 +82,16 @@ for mime in (MIME"text/plain", MIME"text/x-tex", MIME"text/x-latex")
     end
 end
 
+function Base.show(io::IO, ::MIME"text/latex", measure::Complex{<:Measurement})
+    print(io, "\$")
+    Base.show(io, "text/x-latex", measure)
+    print(io, "\$")
+end
+
+function Base.show(io::IO, measure::Complex{<:Measurement})
+    Base.show(io, "text/plain", measure)
+end
+
 # This is adapted from base/show.jl for Complex type.
 function Base.alignment(io::IO, measure::Measurement)
     m = match(r"^(.*[\±])(.*)$", sprint(show, measure, context=io, sizehint=0))
diff --git a/src/utils.jl b/src/utils.jl
index 665b98ff..d1951282 100644
--- a/src/utils.jl
+++ b/src/utils.jl
@@ -20,7 +20,7 @@ computed as the standard score between their difference and 0:
 
     stdscore(measure_1 - measure_2, 0)
 """
-stdscore(a::Measurement{S}, b::Measurement{T}) where {S<:AbstractFloat,T<:AbstractFloat} =
+stdscore(a::Measurement{S}, b::Measurement{T}) where {S<:Real,T<:Real} =
     stdscore(a - b, zero(promote_type(S, T)))
 
 # Weighted Average with Inverse-Variance Weighting
@@ -39,7 +39,7 @@ end
 
 # Derivative and Gradient
 derivative(a::Measurement{F},
-           tag::Tuple{T, T, UInt64}) where {F<:AbstractFloat, T<:AbstractFloat} =
+           tag::Tuple{T, T, UInt64}) where {F<:Real, T<:Real} =
                get(a.der, tag, zero(F))
 
 """
@@ -98,7 +98,7 @@ of an independent `Measurement`, and the value is the absolute value of the
 product between its uncertainty and the partial derivative of `x` with respect
 to this `Measurement`.
 """
-function uncertainty_components(x::Measurement{T}) where {T<:AbstractFloat}
+function uncertainty_components(x::Measurement{T}) where {T<:Real}
     out = Dict{Tuple{T, T, UInt64}, T}()
     for var in keys(x.der)
         out[var] = abs(var[2] * Measurements.derivative(x, var))
diff --git a/test/runtests.jl b/test/runtests.jl
index b5f49867..14f13f83 100644
--- a/test/runtests.jl
+++ b/test/runtests.jl
@@ -7,7 +7,7 @@ else
     Aqua.test_all(Measurements; stale_deps=false)
     Aqua.test_stale_deps(Measurements; ignore=[:RecipesBase, :Requires])
 end
-if VERSION ≥ v"1.6"
+if VERSION ≥ v"1.7"
     Aqua.test_project_toml_formatting(Measurements)
 end
 
@@ -24,8 +24,8 @@ end
 isapprox(x::Measurement, y::Measurement; rest...) =
     isapprox(x.val, y.val; nans = true, rest...) &&
     isapprox(x.err, y.err; nans = true, rest...)
-isapprox(x::Complex{Measurement{<:AbstractFloat}},
-         y::Complex{Measurement{<:AbstractFloat}}; rest...) =
+isapprox(x::Complex{Measurement{<:Real}},
+         y::Complex{Measurement{<:Real}}; rest...) =
              isapprox(real(x), real(y); nans = true, rest...) &&
              isapprox(imag(x), imag(y); nans = true, rest...)
 # This is bit strict, but the idea is that in the tests we want `Measurement`s to be
@@ -784,6 +784,11 @@ end
                        repr("text/plain", [w 10x; 100y 1000w]))
     end
     @testset "Complex numbers" begin
+        @test repr(zz, context=:compact=>true) == "(-0.534±0.031)-(0.534±0.031)im"
+        @test repr(zz) == "(-0.534 ± 0.031) - (0.534 ± 0.031)im"
+        @test repr(complex(x, w)) == "(3.0 ± 0.1) - (0.5 ± 0.03)im"
+        @test repr(complex(w, y)) == "(-0.5 ± 0.03) + (4.0 ± 0.2)im"
+
         @test repr("text/plain", zz, context=:compact=>true) == "(-0.534±0.031)-(0.534±0.031)im"
         @test repr("text/plain", zz) == "(-0.534 ± 0.031) - (0.534 ± 0.031)im"
         @test repr("text/plain", complex(x, w)) == "(3.0 ± 0.1) - (0.5 ± 0.03)im"
@@ -1073,3 +1078,9 @@ end
         @test base_numeric_type(typeof(x)) == T
     end
 end
+
+@static if Base.VERSION >= v"1.9" # We only support Symbolics w/ Julia 1.9+.
+  @testset "Symbolics" begin
+    include("symbolics.jl")
+  end
+end
diff --git a/test/symbolics.jl b/test/symbolics.jl
new file mode 100644
index 00000000..648dc716
--- /dev/null
+++ b/test/symbolics.jl
@@ -0,0 +1,131 @@
+using Symbolics
+
+isapprox(x::Symbolics.Num, y::Symbolics.Num; rest...) = isequal(x, y)
+
+@testset "measurement" begin
+    @variables x, x_err
+
+    @test typeof(@inferred(measurement(x))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, 0))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, 1))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, 0.0))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, 1.0))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, big(0)))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, big(1)))) == Measurement{Num}
+    @test typeof(@inferred(measurement(x, x_err))) == Measurement{Num}
+    @test typeof(@inferred(measurement(0, x_err))) == Measurement{Num}
+    @test typeof(@inferred(measurement(0.0, x_err))) == Measurement{Num}
+    @test typeof(@inferred(measurement(big(0), x_err))) == Measurement{Num}
+
+    @test iszero(@inferred(measurement(x)).err)
+    @test @inferred(measurement(x)).tag === UInt64(0)
+    @test length(@inferred(measurement(x)).der) == 0
+    @test !(@inferred(measurement(x, x_err)).tag === UInt64(0))
+    @test length(@inferred(measurement(x, 0)).der) == 0
+    @test length(@inferred(measurement(x, x_err)).der) == 1
+end
+
+@testset "Type representation" begin
+    @variables x_val, x_err
+
+    x = measurement(x_val, x_err)
+
+    # test pretty printing at the REPL
+    @test_throws ErrorException repr(x, context=:error_digits=>-4)
+    @test repr(x) == "x_val ± x_err"
+    @test repr(x, context=:compact=>true) == "x_val±x_err"
+    for error_digits in (0, 7)
+        @test repr(x, context=:error_digits=>error_digits) == "x_val ± x_err"
+        @test repr(x, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "x_val±x_err"
+    end
+
+    @test repr("text/plain", x) == "x_val ± x_err"
+    @test repr("text/plain", x, context=:compact=>true) == "x_val±x_err"
+    for error_digits in (0, 7)
+        @test repr("text/plain", x, context=:error_digits=>error_digits) == "x_val ± x_err"
+        @test repr("text/plain", x, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "x_val±x_err"
+    end
+
+    @test repr("text/latex", x) == "\$x_val \\pm x_err\$"
+    @test repr("text/latex", x, context=:compact=>true) == "\$x_val\\pmx_err\$"
+    for error_digits in (0, 7)
+        @test repr("text/latex", x, context=:error_digits=>error_digits) == "\$x_val \\pm x_err\$"
+        @test repr("text/latex", x, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "\$x_val\\pmx_err\$"
+    end
+
+    for mime in ("text/x-tex", "text/x-latex")
+        @test repr(mime, x) == "x_val \\pm x_err"
+        @test repr(mime, x, context=:compact=>true) == "x_val\\pmx_err"
+        for error_digits in (0, 7)
+          @test repr(mime, x, context=:error_digits=>error_digits) == "x_val \\pm x_err"
+          @test repr(mime, x, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "x_val\\pmx_err"
+        end
+    end
+
+    @variables y_val, y_err
+    y = measurement(y_val, y_err)
+    c = complex(x, y)
+
+    @test_throws ErrorException repr(c, context=:error_digits=>-4)
+    @test repr(c) == "(x_val ± x_err) + (y_val ± y_err)im"
+    @test repr(c, context=:compact=>true) == "(x_val±x_err)+(y_val±y_err)im"
+    for error_digits in (0, 7)
+        @test repr(c, context=:error_digits=>error_digits) == "(x_val ± x_err) + (y_val ± y_err)im"
+        @test repr(c, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "(x_val±x_err)+(y_val±y_err)im"
+    end
+
+    @test repr("text/plain", c) == "(x_val ± x_err) + (y_val ± y_err)im"
+    @test repr("text/plain", c, context=:compact=>true) == "(x_val±x_err)+(y_val±y_err)im"
+    for error_digits in (0, 7)
+        @test repr("text/plain", c, context=:error_digits=>error_digits) == "(x_val ± x_err) + (y_val ± y_err)im"
+        @test repr("text/plain", c, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "(x_val±x_err)+(y_val±y_err)im"
+    end
+
+    @test repr("text/latex", c) == "\$(x_val \\pm x_err) + (y_val \\pm y_err)im\$"
+    @test repr("text/latex", c, context=:compact=>true) == "\$(x_val\\pmx_err)+(y_val\\pmy_err)im\$"
+    for error_digits in (0, 7)
+        @test repr("text/latex", c, context=:error_digits=>error_digits) == "\$(x_val \\pm x_err) + (y_val \\pm y_err)im\$"
+        @test repr("text/latex", c, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "\$(x_val\\pmx_err)+(y_val\\pmy_err)im\$"
+    end
+
+    for mime in ("text/x-tex", "text/x-latex")
+        @test repr(mime, c) == "(x_val \\pm x_err) + (y_val \\pm y_err)im"
+        @test repr(mime, c, context=:compact=>true) == "(x_val\\pmx_err)+(y_val\\pmy_err)im"
+        for error_digits in (0, 7)
+          @test repr(mime, c, context=:error_digits=>error_digits) == "(x_val \\pm x_err) + (y_val \\pm y_err)im"
+          @test repr(mime, c, context=IOContext(stdout, :error_digits=>error_digits, :compact=>true)) == "(x_val\\pmx_err)+(y_val\\pmy_err)im"
+        end
+    end
+end
+
+##### Mathematical Operations
+@testset "Sin" begin
+    @variables x_val, x_err
+
+    x = measurement(x_val, x_err)
+
+    @test @inferred(sin(x)) ≈ measurement(sin(x_val), abs(x_err*cos(x_val)))
+end
+@testset "Addition" begin
+    @variables x_val, x_err, y_val, y_err
+
+    x = measurement(x_val, x_err)
+    y = measurement(y_val, y_err)
+
+    # abs2(v) === abs(x)^2, but maybe faster.
+    @test @inferred(x + y) ≈ measurement(x_val + y_val, sqrt(abs2(x_err) + abs2(y_err)))
+end
+@testset "Addition w/ correlated arguments" begin
+    @variables x_val, x_err, y_val, y_err
+
+    x = measurement(x_val, x_err)
+    y = measurement(y_val, y_err)
+
+    @variables base_val, base_err
+    base = measurement(base_val, base_err)
+    x += base
+    y += base
+
+    # abs2(v) === abs(x)^2, but maybe faster.
+    @test @inferred(x + y) ≈ measurement(x_val + y_val + 2*base_val, sqrt(abs2(x_err) + abs2(y_err) + abs2(2*base_err)))
+end