diff --git a/src/samplers/gamma.jl b/src/samplers/gamma.jl
index 9a40da98a..1cd62f22b 100644
--- a/src/samplers/gamma.jl
+++ b/src/samplers/gamma.jl
@@ -225,6 +225,6 @@ end
 
 function rand(rng::AbstractRNG, s::GammaIPSampler)
     x = rand(rng, s.s)
-    e = randexp(rng)
+    e = randexp(rng, typeof(x))
     x*exp(s.nia*e)
 end
diff --git a/src/univariate/continuous/exponential.jl b/src/univariate/continuous/exponential.jl
index d1c487057..573a469b0 100644
--- a/src/univariate/continuous/exponential.jl
+++ b/src/univariate/continuous/exponential.jl
@@ -105,7 +105,7 @@ cf(d::Exponential, t::Real) = 1/(1 - t * im * scale(d))
 
 
 #### Sampling
-rand(rng::AbstractRNG, d::Exponential) = xval(d, randexp(rng))
+rand(rng::AbstractRNG, d::Exponential{T}) where {T} = xval(d, randexp(rng, float(T)))
 
 
 #### Fit model
diff --git a/src/univariate/continuous/tdist.jl b/src/univariate/continuous/tdist.jl
index 4e55e8cfb..873925c3e 100644
--- a/src/univariate/continuous/tdist.jl
+++ b/src/univariate/continuous/tdist.jl
@@ -82,7 +82,7 @@ end
 function rand(rng::AbstractRNG, d::TDist)
     ν = d.ν
     z = sqrt(rand(rng, Chisq{typeof(ν)}(ν)) / ν)
-    return randn(rng) / (isinf(ν) ? one(z) : z)
+    return randn(rng, typeof(z)) / (isinf(ν) ? one(z) : z)
 end
 
 function cf(d::TDist{T}, t::Real) where T <: Real
diff --git a/test/univariate/continuous/chisq.jl b/test/univariate/continuous/chisq.jl
index 156bcdb5f..534b7aaf0 100644
--- a/test/univariate/continuous/chisq.jl
+++ b/test/univariate/continuous/chisq.jl
@@ -1,2 +1,9 @@
-test_cgf(Chisq(1), (0.49, -1, -100, -1f6))
-test_cgf(Chisq(3), (0.49, -1, -100, -1f6))
+
+@testset "Chisq" begin
+    test_cgf(Chisq(1), (0.49, -1, -100, -1.0f6))
+    test_cgf(Chisq(3), (0.49, -1, -100, -1.0f6))
+
+    for T in (Float32, Float64)
+        @test @inferred(rand(Chisq(T(1)))) isa T
+    end
+end
diff --git a/test/univariate/continuous/exponential.jl b/test/univariate/continuous/exponential.jl
index 6704554a7..191528fd7 100644
--- a/test/univariate/continuous/exponential.jl
+++ b/test/univariate/continuous/exponential.jl
@@ -1,4 +1,10 @@
 
-test_cgf(Exponential(1), (0.9, -1, -100f0, -1e6))
-test_cgf(Exponential(0.91), (0.9, -1, -100f0, -1e6))
-test_cgf(Exponential(10  ), (0.08, -1, -100f0, -1e6))
+@testset "Exponential" begin
+    test_cgf(Exponential(1), (0.9, -1, -100f0, -1e6))
+    test_cgf(Exponential(0.91), (0.9, -1, -100f0, -1e6))
+    test_cgf(Exponential(10  ), (0.08, -1, -100f0, -1e6))
+
+    for T in (Float32, Float64)
+        @test @inferred(rand(Exponential(T(1)))) isa T
+    end
+end
diff --git a/test/univariate/continuous/gamma.jl b/test/univariate/continuous/gamma.jl
index 3937af191..8c1887786 100644
--- a/test/univariate/continuous/gamma.jl
+++ b/test/univariate/continuous/gamma.jl
@@ -1,25 +1,33 @@
 using Test, Distributions, OffsetArrays
 
-test_cgf(Gamma(1  ,1  ), (0.9, -1, -100f0, -1e6))
-test_cgf(Gamma(10 ,1  ), (0.9, -1, -100f0, -1e6))
-test_cgf(Gamma(0.2, 10), (0.08, -1, -100f0, -1e6))
+@testset "Gamma" begin
+    test_cgf(Gamma(1, 1), (0.9, -1, -100.0f0, -1e6))
+    test_cgf(Gamma(10, 1), (0.9, -1, -100.0f0, -1e6))
+    test_cgf(Gamma(0.2, 10), (0.08, -1, -100.0f0, -1e6))
 
-@testset "Gamma suffstats and OffsetArrays" begin
-    a = rand(Gamma(), 11)
-    wa = 1.0:11.0
+    @testset "Gamma suffstats and OffsetArrays" begin
+        a = rand(Gamma(), 11)
+        wa = 1.0:11.0
 
-    resulta = @inferred(suffstats(Gamma, a))
+        resulta = @inferred(suffstats(Gamma, a))
 
-    resultwa = @inferred(suffstats(Gamma, a, wa))
+        resultwa = @inferred(suffstats(Gamma, a, wa))
 
-    b = OffsetArray(a, -5:5)
-    wb = OffsetArray(wa, -5:5)
+        b = OffsetArray(a, -5:5)
+        wb = OffsetArray(wa, -5:5)
 
-    resultb = @inferred(suffstats(Gamma, b))
-    @test resulta == resultb
+        resultb = @inferred(suffstats(Gamma, b))
+        @test resulta == resultb
 
-    resultwb = @inferred(suffstats(Gamma, b, wb))
-    @test resultwa == resultwb
+        resultwb = @inferred(suffstats(Gamma, b, wb))
+        @test resultwa == resultwb
 
-    @test_throws DimensionMismatch suffstats(Gamma, a, wb)
+        @test_throws DimensionMismatch suffstats(Gamma, a, wb)
+    end
+
+    for T in (Float32, Float64)
+        @test @inferred(rand(Gamma(T(1), T(1)))) isa T
+        @test @inferred(rand(Gamma(1/T(2), T(1)))) isa T
+        @test @inferred(rand(Gamma(T(2), T(1)))) isa T
+    end
 end
diff --git a/test/univariate/continuous/tdist.jl b/test/univariate/continuous/tdist.jl
index 16fab2812..127b99243 100644
--- a/test/univariate/continuous/tdist.jl
+++ b/test/univariate/continuous/tdist.jl
@@ -3,10 +3,17 @@ using ForwardDiff
 
 using Test
 
-@testset "Type stability of `rand` (#1614)" begin
-    if VERSION >= v"1.9.0-DEV.348"
-        # randn(::BigFloat) was only added in https://github.com/JuliaLang/julia/pull/44714
-        @inferred(rand(TDist(big"1.0")))
+@testset "TDist" begin
+    @testset "Type stability of `rand` (#1614)" begin
+        if VERSION >= v"1.9.0-DEV.348"
+            # randn(::BigFloat) was only added in https://github.com/JuliaLang/julia/pull/44714
+            @inferred(rand(TDist(big"1.0")))
+        end
+        @inferred(rand(TDist(ForwardDiff.Dual(1.0))))
+
+    end
+
+    for T in (Float32, Float64)
+        @test @inferred(rand(TDist(T(1)))) isa T
     end
-    @inferred(rand(TDist(ForwardDiff.Dual(1.0))))
 end