title: Understanding how pseudo random numbers in Matlab imply statistical independence
body: <p>Consider the following Matlab code in which I generate some data using pseudo random number generator.
  I would like your help to understand "how" random are these numbers from a statistical point of view, in the terms I explain below.</p>
  
  <p>I first set some parameters</p>
  
  <pre><code>%%%%%%%%Parameters
  clear
  rng default
  Xsup=-1:6; 
  Zsup=1:10; 
  n_m=200; 
  n_w=200; 
  R=n_m;
  </code></pre>
  
  <p>Then I generate the data</p>
  
  <pre><code>%%%%%%%%Creation of data [XZ,etapair,zetapair,etasingle,zetasingle]
  
  %Vector X of dimension n_mx1
  idX=randi(size(Xsup,2),n_m,1); %n_mx1
  X=Xsup(idX).'; %n_mx1
  
  %Vector Z of dimension n_wx1
  idZ=randi(size(Zsup,2),n_w,1); 
  Z=Zsup(idZ).'; %n_wx1
  
  %Combine X and Z in a matrix XZ of dimension (n_m*n_w)x2 
  which lists all possible combinations of values in X and Z
  [cX, cZ] = ndgrid(X,Z);
  XZ = [cX(:), cZ(:)]; %(n_m*n_w)x2
  
  %Vector etapair of dimension (n_m*n_w)x1
  etapair=randn(n_m*n_w,1); %(n_m*n_w)x1
  
  %Vector zetapair of dimension (n_m*n_w)x1
  zetapair=randn(n_m*n_w,1); %(n_m*n_w)x1
  
  %Vector etasingle of dimension (n_m*n_w)x1
  etasingle=max(randn(n_m,R),[],2); %n_mx1 
  etasingle=repmat(etasingle, n_w,1); %(n_m*n_w)x1
  
  %Vector zetasingle of dimension (n_m*n_w)x1
  zetasingle=max(randn(n_w,R),[],2); %n_wx1
  zetasingle=kron(zetasingle, ones(n_m,1)); %(n_m*n_w)x1
  </code></pre>
  
  <p>Let me now translate these draws into statistical terms:</p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>X(t)</code> can be thought as a realisation of a random variable <code>X_t</code></p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>Z(t)</code> can be thought as a realisation of a random variable <code>Z_t</code></p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>etapair(t)</code> can be thought as a realisation of a random variable <code>E_t</code></p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>zetapair(t)</code> can be thought as a realisation of a random variable <code>Q_t</code></p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>etasingle(t)</code> can be thought as a realisation of a random variable <code>Y_t</code></p>
  
  <p>For <code>t=1,...,n_w*n_m</code>, <code>zetasingle(t)</code> can be thought as a realisation of a random variable <code>S_t</code></p>
  
  <p>My belief was that the pseudo random number generator in Matlab allows to claim that
  <code>(X_1,X_2,..., Z_1,Z_2,...,E_1,E_2,..., Q_1,Q_2...,Y_1,Y_2,...,S_1,S_2,...)</code> are mutually independent
  as explained <a href="https://stackoverflow.com/questions/48888423/independence-of-random-variables-in-a-matlab-simulation">here</a></p>
  
  <p>As a check of this hypothetical claim, I define <code>W_t:=-E_t-Q_t+Y_t+S_t</code> and empirically compute <code>Pr(W_t&lt;=1|X_t=5, Z_t=1)</code></p>
  
  <p>If mutual independence holds, then <code>Pr(W_t&lt;=1|X_t=5, Z_t=1)=Pr(W_t&lt;=1)</code> and their empirical counterparts below, named <code>option1</code> and <code>option2</code>, should be ALMOST the same.</p>
  
  <pre><code>%option 1
  num1=zeros(n_m*n_w,1);
  for h=1:n_m*n_w
      if -etapair(h)-zetapair(h)+etasingle(h)+zetasingle(h)&lt;=1 &amp;&amp; XZ(h,1)==5 &amp;&amp; XZ(h,2)==1
          num1(h)=1;
      end
  end
  den1=zeros(n_m*n_w,1);
  for h=1:n_m*n_w
      if  XZ(h,1)==5 &amp;&amp; XZ(h,2)==1
          den1(h)=1;
      end
  end
  option1=sum(num1)/sum(den1);
  
  %option 2
  num2=zeros(n_m*n_w,1);
  for h=1:n_m*n_w
      if -etapair(h)-zetapair(h)+etasingle(h)+zetasingle(h)&lt;=1 
          num2(h)=1;
      end
  end
  option2=sum(num2)/(n_m*n_w);
  </code></pre>
  
  <p><strong>Question</strong>: the difference between <code>option1</code> (=0.0021) and <code>option2</code> (=0.0012) is referred to the "ALMOST" or I am doing something wrong?</p>