diff --git a/src/Numerics/Generate.cs b/src/Numerics/Generate.cs
index 69f11ef3..b447e3ae 100644
--- a/src/Numerics/Generate.cs
+++ b/src/Numerics/Generate.cs
@@ -589,7 +589,7 @@ namespace MathNet.Numerics
         /// </summary>
         public static double[] WhiteGaussianNoise(int length, double mean, double standardDeviation)
         {
-            return Normal.Samples(MersenneTwister.Default, mean, standardDeviation).Take(length).ToArray();
+            return Normal.Samples(SystemRandomSource.Default, mean, standardDeviation).Take(length).ToArray();
         }
 
         /// <summary>
@@ -597,7 +597,7 @@ namespace MathNet.Numerics
         /// </summary>
         public static IEnumerable<double> WhiteGaussianNoiseSequence(double mean, double standardDeviation)
         {
-            return Normal.Samples(MersenneTwister.Default, mean, standardDeviation);
+            return Normal.Samples(SystemRandomSource.Default, mean, standardDeviation);
         }
 
         /// <summary>
@@ -610,7 +610,7 @@ namespace MathNet.Numerics
         /// <param name="location">Location mu-parameter of the stable distribution</param>
         public static double[] StableNoise(int length, double alpha, double beta, double scale, double location)
         {
-            return Stable.Samples(MersenneTwister.Default, alpha, beta, scale, location).Take(length).ToArray();
+            return Stable.Samples(SystemRandomSource.Default, alpha, beta, scale, location).Take(length).ToArray();
         }
 
         /// <summary>
@@ -622,7 +622,7 @@ namespace MathNet.Numerics
         /// <param name="location">Location mu-parameter of the stable distribution</param>
         public static IEnumerable<double> StableNoiseSequence(double alpha, double beta, double scale, double location)
         {
-            return Stable.Samples(MersenneTwister.Default, alpha, beta, scale, location);
+            return Stable.Samples(SystemRandomSource.Default, alpha, beta, scale, location);
         }
     }
 }
diff --git a/src/Numerics/LinearAlgebra/Builder.cs b/src/Numerics/LinearAlgebra/Builder.cs
index 3cc23640..9d5251bf 100644
--- a/src/Numerics/LinearAlgebra/Builder.cs
+++ b/src/Numerics/LinearAlgebra/Builder.cs
@@ -34,7 +34,6 @@ using System.Linq;
 using MathNet.Numerics.Distributions;
 using MathNet.Numerics.LinearAlgebra.Solvers;
 using MathNet.Numerics.LinearAlgebra.Storage;
-// TODO: split up and move to proper folders
 using MathNet.Numerics.Random;
 
 namespace MathNet.Numerics.LinearAlgebra.Double
@@ -350,7 +349,6 @@ namespace MathNet.Numerics.LinearAlgebra
     using Complex64 = Numerics.Complex;
 #else
     using Complex64 = System.Numerics.Complex;
-
 #endif
 
     /// <summary>
@@ -447,7 +445,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// </summary>
         public Matrix<T> Random(int rows, int columns)
         {
-            return Random(rows, columns, new Normal(MersenneTwister.Default));
+            return Random(rows, columns, new Normal(SystemRandomSource.Default));
         }
 
         /// <summary>
@@ -455,7 +453,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// </summary>
         public Matrix<T> Random(int rows, int columns, int seed)
         {
-            return Random(rows, columns, new Normal(new MersenneTwister(seed, true)));
+            return Random(rows, columns, new Normal(new SystemRandomSource(seed, true)));
         }
 
         /// <summary>
@@ -474,7 +472,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// </summary>
         public Matrix<T> RandomPositiveDefinite(int order)
         {
-            var a = Random(order, order, new Normal(MersenneTwister.Default));
+            var a = Random(order, order, new Normal(SystemRandomSource.Default));
             return a.ConjugateTransposeThisAndMultiply(a);
         }
 
@@ -484,7 +482,7 @@ namespace MathNet.Numerics.LinearAlgebra
         /// </summary>
         public Matrix<T> RandomPositiveDefinite(int order, int seed)
         {
-            var a = Random(order, order, new Normal(new MersenneTwister(seed, true)));
+            var a = Random(order, order, new Normal(new SystemRandomSource(seed, true)));
             return a.ConjugateTransposeThisAndMultiply(a);
         }
 
@@ -1395,15 +1393,15 @@ namespace MathNet.Numerics.LinearAlgebra
         /// </summary>
         public Vector<T> Random(int length)
         {
-            return Random(length, new Normal(MersenneTwister.Default));
+            return Random(length, new Normal(SystemRandomSource.Default));
         }
 
         /// <summary>
-        /// Create a new dense vector with values sampled from the standard distribution with a mersenne twister random source.
+        /// Create a new dense vector with values sampled from the standard distribution with a system random source.
         /// </summary>
         public Vector<T> Random(int length, int seed)
         {
-            return Random(length, new Normal(new MersenneTwister(seed, true)));
+            return Random(length, new Normal(new SystemRandomSource(seed, true)));
         }
 
         /// <summary>
diff --git a/src/Numerics/Statistics/MCMC/HybridMC.cs b/src/Numerics/Statistics/MCMC/HybridMC.cs
index 61411197..5457edee 100644
--- a/src/Numerics/Statistics/MCMC/HybridMC.cs
+++ b/src/Numerics/Statistics/MCMC/HybridMC.cs
@@ -85,7 +85,7 @@ namespace MathNet.Numerics.Statistics.Mcmc
         /// <param name="burnInterval">The number of iterations in between returning samples.</param>
         /// <exception cref="ArgumentOutOfRangeException">When the number of burnInterval iteration is negative.</exception>
         public HybridMC(double[] x0, DensityLn<double[]> pdfLnP, int frogLeapSteps, double stepSize, int burnInterval = 0)
-            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, new double[x0.Count()], MersenneTwister.Default, Grad)
+            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, new double[x0.Count()], SystemRandomSource.Default, Grad)
         {
             for (int i = 0; i < _length; i++)
             {
@@ -109,7 +109,7 @@ namespace MathNet.Numerics.Statistics.Mcmc
         /// the components of the momentum.</param>
         /// <exception cref="ArgumentOutOfRangeException">When the number of burnInterval iteration is negative.</exception>
         public HybridMC(double[] x0, DensityLn<double[]> pdfLnP, int frogLeapSteps, double stepSize, int burnInterval, double[] pSdv)
-            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, pSdv, MersenneTwister.Default)
+            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, pSdv, SystemRandomSource.Default)
         {
         }
 
diff --git a/src/Numerics/Statistics/MCMC/MCMCSampler.cs b/src/Numerics/Statistics/MCMC/MCMCSampler.cs
index a23f3411..017a0ac8 100644
--- a/src/Numerics/Statistics/MCMC/MCMCSampler.cs
+++ b/src/Numerics/Statistics/MCMC/MCMCSampler.cs
@@ -107,7 +107,7 @@ namespace MathNet.Numerics.Statistics.Mcmc
         {
             Accepts = 0;
             Samples = 0;
-            RandomSource = MersenneTwister.Default;
+            RandomSource = SystemRandomSource.Default;
         }
 
         /// <summary>
diff --git a/src/Numerics/Statistics/MCMC/UnivariateHybridMC.cs b/src/Numerics/Statistics/MCMC/UnivariateHybridMC.cs
index 0c1f9f8a..1cd72c24 100644
--- a/src/Numerics/Statistics/MCMC/UnivariateHybridMC.cs
+++ b/src/Numerics/Statistics/MCMC/UnivariateHybridMC.cs
@@ -83,7 +83,7 @@ namespace MathNet.Numerics.Statistics.Mcmc
         /// the momentum.</param>
         /// <exception cref="ArgumentOutOfRangeException">When the number of burnInterval iteration is negative.</exception>
         public UnivariateHybridMC(double x0, DensityLn<double> pdfLnP, int frogLeapSteps, double stepSize, int burnInterval = 0, double pSdv = 1)
-            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, pSdv, MersenneTwister.Default)
+            : this(x0, pdfLnP, frogLeapSteps, stepSize, burnInterval, pSdv, SystemRandomSource.Default)
         {
         }
 
diff --git a/src/UnitTests/DistributionTests/CommonDistributionTests.cs b/src/UnitTests/DistributionTests/CommonDistributionTests.cs
index 8fcb6862..843c4314 100644
--- a/src/UnitTests/DistributionTests/CommonDistributionTests.cs
+++ b/src/UnitTests/DistributionTests/CommonDistributionTests.cs
@@ -183,7 +183,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
         [Test]
         public void SampleFollowsCorrectDistribution()
         {
-            Random rnd = new MersenneTwister(1);
+            Random rnd = new SystemRandomSource(1);
 
             foreach (var dd in _discreteDistributions)
             {
@@ -216,7 +216,7 @@ namespace MathNet.Numerics.UnitTests.DistributionTests
         [Test]
         public void SamplesFollowsCorrectDistribution()
         {
-            Random rnd = new MersenneTwister(1);
+            Random rnd = new SystemRandomSource(1);
 
             foreach (var dd in _discreteDistributions)
             {
diff --git a/src/UnitTests/InterpolationTests/LinearInterpolationCase.cs b/src/UnitTests/InterpolationTests/LinearInterpolationCase.cs
index 4e989c90..1a7a5da9 100644
--- a/src/UnitTests/InterpolationTests/LinearInterpolationCase.cs
+++ b/src/UnitTests/InterpolationTests/LinearInterpolationCase.cs
@@ -52,7 +52,7 @@ namespace MathNet.Numerics.UnitTests.InterpolationTests
         public static void Build(out double[] x, out double[] y, out double[] xtest, out double[] ytest, int samples = 3, double sampleOffset = -0.5, double slope = 2.0, double intercept = -1.0)
         {
             // Fixed-seed "random" distribution to ensure we always test with the same data
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister(42));
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(42));
 
             // build linear samples
             x = new double[samples];
diff --git a/src/UnitTests/LinearAlgebraTests/Complex/DiagonalMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex/DiagonalMatrixTests.cs
index 4c9b547e..2ed3c9a5 100644
--- a/src/UnitTests/LinearAlgebraTests/Complex/DiagonalMatrixTests.cs
+++ b/src/UnitTests/LinearAlgebraTests/Complex/DiagonalMatrixTests.cs
@@ -368,7 +368,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DenseDiagonalMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex>.Build.Random(8, 3, dist);
@@ -385,7 +385,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DenseDiagonalMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex>.Build.Random(8, 3, dist);
@@ -402,7 +402,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DenseDiagonalMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex>.Build.Random(3, 8, dist);
@@ -419,7 +419,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DiagonalDenseMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex>.Build.Random(3, 8, dist);
@@ -436,7 +436,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DiagonalDenseMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex>.Build.Random(8, 3, dist);
@@ -453,7 +453,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex
         [Test]
         public void DiagonalDenseMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex>.Build.Random(3, 8, dist);
diff --git a/src/UnitTests/LinearAlgebraTests/Complex32/DiagonalMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Complex32/DiagonalMatrixTests.cs
index 7d83586e..89589d70 100644
--- a/src/UnitTests/LinearAlgebraTests/Complex32/DiagonalMatrixTests.cs
+++ b/src/UnitTests/LinearAlgebraTests/Complex32/DiagonalMatrixTests.cs
@@ -364,7 +364,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DenseDiagonalMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex32>.Build.Random(8, 3, dist);
@@ -381,7 +381,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DenseDiagonalMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex32>.Build.Random(8, 3, dist);
@@ -398,7 +398,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DenseDiagonalMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex32>.Build.Random(3, 8, dist);
@@ -415,7 +415,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DiagonalDenseMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex32>.Build.Random(3, 8, dist);
@@ -432,7 +432,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DiagonalDenseMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<Complex32>.Build.Random(8, 3, dist);
@@ -449,7 +449,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Complex32
         [Test]
         public void DiagonalDenseMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<Complex32>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<Complex32>.Build.Random(3, 8, dist);
diff --git a/src/UnitTests/LinearAlgebraTests/Double/DiagonalMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Double/DiagonalMatrixTests.cs
index e1febf20..829f55a2 100644
--- a/src/UnitTests/LinearAlgebraTests/Double/DiagonalMatrixTests.cs
+++ b/src/UnitTests/LinearAlgebraTests/Double/DiagonalMatrixTests.cs
@@ -396,7 +396,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DenseDiagonalMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<double>.Build.Random(8, 3, dist);
@@ -413,7 +413,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DenseDiagonalMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<double>.Build.Random(8, 3, dist);
@@ -430,7 +430,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DenseDiagonalMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<double>.Build.Random(3, 8, dist);
@@ -447,7 +447,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DiagonalDenseMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<double>.Build.Random(3, 8, dist);
@@ -464,7 +464,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DiagonalDenseMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<double>.Build.Random(8, 3, dist);
@@ -481,7 +481,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Double
         [Test]
         public void DiagonalDenseMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<double>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<double>.Build.Random(3, 8, dist);
diff --git a/src/UnitTests/LinearAlgebraTests/Single/DiagonalMatrixTests.cs b/src/UnitTests/LinearAlgebraTests/Single/DiagonalMatrixTests.cs
index 36c7e15a..eb651e9a 100644
--- a/src/UnitTests/LinearAlgebraTests/Single/DiagonalMatrixTests.cs
+++ b/src/UnitTests/LinearAlgebraTests/Single/DiagonalMatrixTests.cs
@@ -362,7 +362,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DenseDiagonalMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<float>.Build.Random(8, 3, dist);
@@ -379,7 +379,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DenseDiagonalMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<float>.Build.Random(8, 3, dist);
@@ -396,7 +396,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DenseDiagonalMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<float>.Build.Random(3, 8, dist);
@@ -413,7 +413,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DiagonalDenseMatrixMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<float>.Build.Random(3, 8, dist);
@@ -430,7 +430,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DiagonalDenseMatrixTransposeAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var tall = Matrix<float>.Build.Random(8, 3, dist);
@@ -447,7 +447,7 @@ namespace MathNet.Numerics.UnitTests.LinearAlgebraTests.Single
         [Test]
         public void DiagonalDenseMatrixTransposeThisAndMultiply()
         {
-            var dist = new ContinuousUniform(-1.0, 1.0, new MersenneTwister());
+            var dist = new ContinuousUniform(-1.0, 1.0, new SystemRandomSource(1));
             Assert.IsInstanceOf<DiagonalMatrix>(Matrix<float>.Build.DiagonalIdentity(3, 3));
 
             var wide = Matrix<float>.Build.Random(3, 8, dist);
diff --git a/src/UnitTests/StatisticsTests/MCMCTests/MetropolisHastingsSamplerTests.cs b/src/UnitTests/StatisticsTests/MCMCTests/MetropolisHastingsSamplerTests.cs
index 90acd16c..896a4126 100644
--- a/src/UnitTests/StatisticsTests/MCMCTests/MetropolisHastingsSamplerTests.cs
+++ b/src/UnitTests/StatisticsTests/MCMCTests/MetropolisHastingsSamplerTests.cs
@@ -49,7 +49,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void MetropolisHastingsConstructor()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisHastingsSampler<double>(0.2, normal.Density, (x, y) => Normal.PDF(x, 0.1, y), x => Normal.Sample(rnd, x, 0.1), 10)
                 {
@@ -68,7 +68,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void SampleTest()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisHastingsSampler<double>(0.2, normal.Density, (x, y) => Normal.PDF(x, 0.1, y), x => Normal.Sample(rnd, x, 0.1), 10)
                 {
@@ -85,7 +85,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void SampleArrayTest()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisHastingsSampler<double>(0.2, normal.Density, (x, y) => Normal.PDF(x, 0.1, y), x => Normal.Sample(rnd, x, 0.1), 10)
                 {
diff --git a/src/UnitTests/StatisticsTests/MCMCTests/MetropolisSamplerTests.cs b/src/UnitTests/StatisticsTests/MCMCTests/MetropolisSamplerTests.cs
index 8c82a726..c6d6756e 100644
--- a/src/UnitTests/StatisticsTests/MCMCTests/MetropolisSamplerTests.cs
+++ b/src/UnitTests/StatisticsTests/MCMCTests/MetropolisSamplerTests.cs
@@ -51,7 +51,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void MetropolisConstructor()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisSampler<double>(0.2, normal.Density, x => Normal.Sample(rnd, x, 0.1), 10);
             Assert.IsNotNull(ms.RandomSource);
@@ -67,7 +67,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void SampleTest()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisSampler<double>(0.2, normal.Density, x => Normal.Sample(rnd, x, 0.1), 10)
                 {
@@ -84,7 +84,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         public void SampleArrayTest()
         {
             var normal = new Normal(0.0, 1.0);
-            var rnd = new MersenneTwister();
+            var rnd = new SystemRandomSource(1);
 
             var ms = new MetropolisSampler<double>(0.2, normal.Density, x => Normal.Sample(rnd, x, 0.1), 10)
                 {
diff --git a/src/UnitTests/StatisticsTests/MCMCTests/RejectionSamplerTests.cs b/src/UnitTests/StatisticsTests/MCMCTests/RejectionSamplerTests.cs
index cdf7a2d0..1b2f8084 100644
--- a/src/UnitTests/StatisticsTests/MCMCTests/RejectionSamplerTests.cs
+++ b/src/UnitTests/StatisticsTests/MCMCTests/RejectionSamplerTests.cs
@@ -48,7 +48,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         [Test]
         public void RejectTest()
         {
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister());
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(1));
             var rs = new RejectionSampler<double>(x => Math.Pow(x, 1.7)*Math.Pow(1.0 - x, 5.3), x => 0.021, uniform.Sample);
             Assert.IsNotNull(rs.RandomSource);
 
@@ -62,7 +62,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         [Test]
         public void SampleTest()
         {
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister());
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(1));
 
             var rs = new RejectionSampler<double>(x => Math.Pow(x, 1.7)*Math.Pow(1.0 - x, 5.3), x => 0.021, uniform.Sample)
                 {
@@ -78,7 +78,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         [Test]
         public void SampleArrayTest()
         {
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister());
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(1));
 
             var rs = new RejectionSampler<double>(x => Math.Pow(x, 1.7)*Math.Pow(1.0 - x, 5.3), x => 0.021, uniform.Sample)
                 {
@@ -94,7 +94,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         [Test]
         public void NoUpperBound()
         {
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister());
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(1));
             var rs = new RejectionSampler<double>(x => Math.Pow(x, 1.7)*Math.Pow(1.0 - x, 5.3), x => Double.NegativeInfinity, uniform.Sample);
             Assert.Throws<ArgumentOutOfRangeException>(() => rs.Sample());
         }
@@ -105,7 +105,7 @@ namespace MathNet.Numerics.UnitTests.StatisticsTests.McmcTests
         [Test]
         public void NullRandomNumberGenerator()
         {
-            var uniform = new ContinuousUniform(0.0, 1.0, new MersenneTwister());
+            var uniform = new ContinuousUniform(0.0, 1.0, new SystemRandomSource(1));
             var rs = new RejectionSampler<double>(x => Math.Pow(x, 1.7)*Math.Pow(1.0 - x, 5.3), x => Double.NegativeInfinity, uniform.Sample);
             Assert.Throws<ArgumentNullException>(() => rs.RandomSource = null);
         }