1

The problem: If a large number of fair N-sided dice are rolled, the average of the simulated rolls is likely to be close to the mean of 1,2,...N i.e. the expected value of one die. For example, the expected value of a 6-sided die is 3.5.

Given N, simulate 1e8 N-sided dice rolls by creating a vector of 1e8 uniformly distributed random integers. Return the difference between the mean of this vector and the mean of integers from 1 to N.


My code:

function dice_diff = loln(N)
       % the mean of integer from 1 to N    
       A = 1:N
       meanN = sum(A)/N;
% I do not have any idea what I am doing here!
       V = randi(1e8);
       meanvector = V/1e8;
       dice_diff = meanvector - meanN;
end
1
  • 5
    What is the question you have? Feb 26, 2020 at 11:48

2 Answers 2

1

First of all, make sure everytime you ask a question that it is as clear as possible, to make it easier for other users to read. If you check how randi works, you can see this:

R = randi(IMAX,N) returns an N-by-N matrix containing pseudorandom integer values drawn from the discrete uniform distribution on 1:IMAX. randi(IMAX,M,N) or randi(IMAX,[M,N]) returns an M-by-N matrix. randi(IMAX,M,N,P,...) or randi(IMAX,[M,N,P,...]) returns an M-by-N-by-P-by-... array. randi(IMAX) returns a scalar. randi(IMAX,SIZE(A)) returns an array the same size as A.

So, if you want to use randi in your problem, you have to use it like this:

V=randi(N, 1e8,1);

and you need some more changes:

function dice_diff = loln(N)
%the mean of integer from 1 to N    
A = 1:N;
meanN = mean(A); 
V = randi(N, 1e8,1);
meanvector = mean(V);
dice_diff = meanvector - meanN;
end

For future problems, try using the command

help randi

And matlab will explain how the function randi (or other function) works.

Make sure to check if the code above gives the desired result

0
1

As pointed out, take a closer look at the use of randi(). From the general case

X = randi([LowerInt,UpperInt],NumRows,NumColumns);  % UpperInt > LowerInt

you can adapt to dice rolling by

Rolls = randi([1 NumSides],NumRolls,NumSamplePaths); 

as an example. Exchanging NumRolls and NumSamplePaths will yield Rolls.', or transpose(Rolls).

According to the Law of Large Numbers, the updated sample average after each roll should converge to the true mean, ExpVal (short for expected value), as the number of rolls (trials) increases. Notice that as NumRolls gets larger, the sample mean converges to the true mean. The image below shows this for two sample paths.

Top Image showing sample mean converging per Law of Large Numbers. Bottom figure showing mean error converging to zero.

To get the sample mean for each number of dice rolls, I used arrayfun() with

CumulativeAvg1 = arrayfun(@(jj)mean(Rolls(1:jj,1)),[1:NumRolls]);

which is equivalent to using the cumulative sum, cumsum(), to get the same result.

CumulativeAvg1 = (cumsum(Rolls(:,1))./(1:NumRolls).');    % equivalent

% MATLAB R2019a
% Create Dice
NumSides = 6;   % positive nonzero integer
NumRolls = 200;
NumSamplePaths = 2;

% Roll Dice
Rolls = randi([1 NumSides],NumRolls,NumSamplePaths);  

% Output Statistics
ExpVal = mean(1:NumSides);
CumulativeAvg1 = arrayfun(@(jj)mean(Rolls(1:jj,1)),[1:NumRolls]);
CumulativeAvgError1 = CumulativeAvg1 - ExpVal;
CumulativeAvg2 = arrayfun(@(jj)mean(Rolls(1:jj,2)),[1:NumRolls]);
CumulativeAvgError2 = CumulativeAvg2 - ExpVal;

% Plot
figure
subplot(2,1,1), hold on, box on
    plot(1:NumRolls,CumulativeAvg1,'b--','LineWidth',1.5,'DisplayName','Sample Path 1')
    plot(1:NumRolls,CumulativeAvg2,'r--','LineWidth',1.5,'DisplayName','Sample Path 2')
    yline(ExpVal,'k-')
    title('Average')
    xlabel('Number of Trials')
    ylim([1 NumSides])
subplot(2,1,2), hold on, box on
    plot(1:NumRolls,CumulativeAvgError1,'b--','LineWidth',1.5,'DisplayName','Sample Path 1')
    plot(1:NumRolls,CumulativeAvgError2,'r--','LineWidth',1.5,'DisplayName','Sample Path 2')
    yline(0,'k-')
    title('Error')
    xlabel('Number of Trials')

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service and acknowledge you have read our privacy policy.

Not the answer you're looking for? Browse other questions tagged or ask your own question.