Perfect random numbers - How can I refine an initial dataset from randn to be perfrectly normal?

2 次查看(过去 30 天)
Hi, I would like to create "perfect" random numbers. By perfect I mean a set of random values that has
a mean of 0.0000, a std of 1.0000, a skewness of 0.0000, and a kurtosis of 3.0000
while each of numbers created with randn comes very close, it is not perfect in that sense. How can I further refine the initial dataset from randn to be "perfectly" normal?
Thanks

回答(5 个)

Teja Muppirala
Teja Muppirala 2013-5-22
If you want to generate a set of numbers that "looks random", has 0 mean, 1 variance, 0 skewness, 3 kurtosis, you can start with an initial dataset using RANDN and then kind of iteratively perturb it until you get what you want.
rng(0)
R = randn(1000,1);
R0 = R;
R_prev= 0;
opts = optimset('TolFun',0,'Display','Off');
smallNumber = 1e-12;
while norm(R-R_prev)/norm(R_prev) > smallNumber;
R_prev = R;
R = R-mean(R); %Fix Mean
R = R/std(R); %Fix Std. Dev
% Fix Skewness
x = fzero( @(x) skewness( R .* exp(x*R) ) ,0 ,opts);
R = R .* exp(x*R);
% Fix Kurtosis
x = fzero( @(x) kurtosis(sign(R) .* (abs(R).^x))-3 ,1 ,opts);
R = sign(R) .* (abs(R).^x);
end
figure, plot(R0,'b');
hold on;
plot(R,'r');
legend({'Original', 'Slightly Modified'});
format long
disp('Original:')
[mean(R0) std(R0) skewness(R0) kurtosis(R0)]
disp('Modified:')
[mean(R) std(R) skewness(R) kurtosis(R)]
  5 个评论
Sven
Sven 2013-6-12
Teja, I have changed the above to be part of a custom formula that reads like follows:
function [vRand ] = fPerfectRandomFigures(iNumber,dStd,dSkew,dKurt)
%UNTITLED Summary of this function goes here
% Detailed explanation goes here
rng('shuffle')
R = randn(iNumber,1);
R_prev= 0;
opts = optimset('TolFun',0,'Display','Off');
smallNumber = 1e-16;
while norm(R-R_prev)/norm(R_prev) > smallNumber;
R_prev = R;
R = R-mean(R); %Fix Mean
R = R/std(R) * dStd; %Fix Std. Dev
x = fzero( @(x) skewness( R .* exp(x*R) )-dSkew ,0 ,opts); % Fix Skewness
R = R .* exp(x*R); % Fix Kurtosis
x = fzero( @(x) kurtosis(sign(R) .* (abs(R).^x))-(3+dKurt) ,1 ,opts);
R = sign(R) .* (abs(R).^x);
end
vRand=R;
end
I now would like to extend this a bit further to also allow to control for the serial correlation of the random numbers (1st order serial correlation). I would like to add a new parameter that can take any value between -1 and 1 and the vector of random numbers would be adjusted to have a serial correlation of that value. How would I do that?
Thanks for any hint...
Sven
süleymand demir
süleymand demir 2020-11-16
Hi Teja,
I wanted to use this code to genarate some datas. But I could not generate data for some conditions.
For example: a mean of 0.0000, a std of 1.0000, a skewness of 1.0000, and a kurtosis of 4.0000.
Could you help me please.

请先登录,再进行评论。


Walter Roberson
Walter Roberson 2013-5-21
You cannot do this without constraining two of the values to be non-random (values that fix the mean and std.)
The Normal Distribution is defined by an infinite process; if you use a finite number of values, you can only approximate the Normal Distribution.
  3 个评论
Walter Roberson
Walter Roberson 2013-5-22
Roger has a FEX contribution for a "random" set of numbers with a fixed sum. That technique can be used with the fixed sum being 0. Call that T. Then T/std(T) will have a sum of 0 and a std of 1.
I have difficulties accepting that the numbers are truly "random" if they have a fixed sum, but if Roger finds it meaningful then I consider that perhaps there is something lacking in my understanding.
Jan
Jan 2013-5-24
编辑:Jan 2013-5-24
I think, Roger applies a specific definition of "random" here: The maximal achievable indepedence between the values might be a solid base.
Considering the total interconnectedness of the properties in the universe, random is only vaguely defined at all. Einstein, Podolsky and Rosen hesitated to accept, that the spins of two conjugated photons can be oriented randomly although the sum of the spins will vanish after the measurement.

请先登录,再进行评论。


Iain
Iain 2013-5-22
It isn't possible.
Just looking at the mean & std, if you do manage to get the mean to be exactly 0 and the std to be exactly 1, how can you possibly add another sample without either changing the mean (and retaining the std), or changing the std (and retaining the mean). Skewness and kurtosis just make it even harder.

Sven
Sven 2013-5-22
To all that commented: thanks for your comments. I guess I should have been a bit more precise in my question regarding being perfect: I was referring to almost perfect only...errors in the 3rd or 4th decimal point as by using this code is absolutely fine for what I am after.
  3 个评论
Sven
Sven 2013-5-24
while it is certainly possible to extend the number of draws to a very high number, I am trying to keep the size of the sample rather small.
Matt J
Matt J 2013-6-15
编辑:Matt J 2013-6-15
Whether the number of draws is large or small is irrelevant. The statistical mean and other moments are what they are independently of how you sample them. My test with large draws was only to measure the statistical moments better and show you that they have the values they should.
If you're saying that you want the empirical moments to deviate from the statistical moments by a limited amount for a certain sample size, then you have no control over that, I'm afraid. The empirical moments are random variables and will have a variance that depends on your sample size and the statistical variance. You cannot change the empirical moments' variances without first changing the statistical distribution of the things you're sampling.

请先登录,再进行评论。


Iain
Iain 2013-6-14
Hmmm:
Take the inverse cdf, at regular intervals.
rnds = icdf('norm',(1/samples):(1/samples):(1-1/samples), 0 , 1);
These numbers are not random, but they, as a set, have a mean of 0, std of close to 1, a skewness of 0, and a kurtosis of "close to" 3.
Generate a random order:
[~, I] = sort(randn(size(rnds)));
randoms = rnds(I);
  2 个评论
Sven
Sven 2013-6-14
Thanks lain, the first part works fine for generating the initial set of random numbers, but I am not sure I follow how the second part of your suggested answer would enforce a prespecified autocorrelation on the vector (see my comment from June 12th above), or was that not part of your answer? Sven

请先登录,再进行评论。

Community Treasure Hunt

Find the treasures in MATLAB Central and discover how the community can help you!

Start Hunting!

Translated by