www.gusucode.com > stats 源码程序 matlab案例代码 > stats/GenerateAQuasiRandomPointSetExample.m

    %% Generate a Quasi-Random Point Set
% This example shows how to use |haltonset| to construct a 2-D Halton 
% quasi-random point set.

% Copyright 2015 The MathWorks, Inc.


%%
% Create a |haltonset| object |p|, that skips the first 1000 values of the
% sequence and then retains every 101st point.
rng default  % For reproducibility
p = haltonset(2,'Skip',1e3,'Leap',1e2)

%%
% The object |p| encapsulates properties of the specified quasi-random
% sequence. The point set is finite, with a length determined by the |Skip|
% and |Leap| properties and by limits on the size of point set indices.

%%
% Use |scramble| to apply reverse-radix scrambling.
p = scramble(p,'RR2')

%%
% Use |net| to generate the first 500 points.
X0 = net(p,500);

%%
% This is equivalent to
X0 = p(1:500,:);

%%
% Values of the point set |X0| are not generated and stored in memory until
% you access |p| using |net| or parenthesis indexing.

%%
% To appreciate the nature of quasi-random numbers, create a scatter plot
% of the two dimensions in |X0|.
scatter(X0(:,1),X0(:,2),5,'r')
axis square
title('{\bf Quasi-Random Scatter}')

%%
% Compare this to a scatter of uniform pseudorandom numbers generated by the
% |rand| function.
X = rand(500,2);
scatter(X(:,1),X(:,2),5,'b')
axis square
title('{\bf Uniform Random Scatter}')

%%
% The quasi-random scatter appears more uniform, avoiding the clumping in
% the pseudorandom scatter.

%%
% In a statistical sense, quasi-random numbers are too uniform to pass
% traditional tests of randomness. For example, a Kolmogorov-Smirnov test,
% performed by |kstest|, is used to assess whether or not a point set has a
% uniform random distribution. When performed repeatedly on uniform
% pseudorandom samples, such as those generated by |rand|, the test produces
% a uniform distribution of _p_-values.
nTests = 1e5;
sampSize = 50;
PVALS = zeros(nTests,1);
for test = 1:nTests
    x = rand(sampSize,1);
    [h,pval] = kstest(x,[x,x]);
    PVALS(test) = pval;
end

histogram(PVALS,100)
h = findobj(gca,'Type','patch');
xlabel('{\it p}-values')
ylabel('Number of Tests')

%%
% The results are quite different when the test is performed repeatedly on
% uniform quasi-random samples.
p = haltonset(1,'Skip',1e3,'Leap',1e2);
p = scramble(p,'RR2');

nTests = 1e5;
sampSize = 50;
PVALS = zeros(nTests,1);
for test = 1:nTests
    x = p(test:test+(sampSize-1),:);
    [h,pval] = kstest(x,[x,x]);
    PVALS(test) = pval;
end

histogram(PVALS,100)
xlabel('{\it p}-values')
ylabel('Number of Tests')

%%
% Small _p_-values call into question the null hypothesis that the data are
% uniformly distributed. If the hypothesis is true, about 5% of the _p_-values
% are expected to fall below 0.05. The results are remarkably consistent in
% their failure to challenge the hypothesis.