Different results when calling a function and writing function in the same script

Hi,

I am very new with matlab so this may be a very basic mistake that I am making. I am running two codes that do the same thing, but they return different result. The first code calls two functions which I wrote in a different script while the second code includes the function in the script itself (ie I do not call a function). The following is the lines of codes from the first code, along with the functions.

% Set up environment
rho = 0.5;
alpha = 0;
T = 1100;
sigma = 1;

% Simulate and run OLS 100 times and save estimate
rho_hat = zeros(100,1);

for i = 1:100 
    
    % Simulate AR and subset data
    Y = sim_AR(rho,alpha,sigma,T);
    Yt_100 = Y(1001:end);
    Ytmin1_100 = Y(1000:end-1);
    Yt_all = [Yt_100 Ytmin1_100];
    
    [beta,cvar,se] = calc_OLS(Yt_100,Ytmin1_100);
    
    rho_hat(i) = beta;
end

mean = mean(rho_hat);

As you can see, these codes call two functions which are the following

function [Y] = sim_AR(alpha,rho,sigma,T)
    y_0 = 0;
    Y = [];
    
    for t = 1:T
        eps = normrnd(0,sigma);
        if t == 1
            y = alpha + rho * y_0 + eps;
        else
            y = alpha + rho * Y(t-1) + eps;
        end
        
        Y = [Y ; y];
        
    end
end

function [betahat, cvar, se] = calc_OLS(Y,X)
    sz = size(X);
    n = sz(1);
    k = sz(2);
    betahat = inv(X.' * X) * X.' * Y;
    epsilon = Y - X * betahat;
    cvar = epsilon.' * epsilon / (n - k);
    var = diag(cvar);
    se = var.^(1/2);
end

The second code that I run is basically putting the codes of the function in one long code

    clear all
    clc

    rho = 0.5;
    alpha = 0;
    T = 1100;
    sigma = 1;
    
    rho_hat = [];

for i = 1:100 
    
    y_0 = 0;
    Y = [];
    for t = 1:T
        eps = normrnd(0,sigma);
        if t == 1
            y = alpha + rho * y_0 + eps;
        else
            y = alpha + rho * Y(t-1) + eps;
        end
        
        Y = [Y ; y];
        
    end
    
    Yt_100 = Y(1001:end);
    Ytmin1_100 = Y(1000:end-1);
    Yt_all = [Yt_100 Ytmin1_100];    
    
        [beta,cvar,se] = calc_OLS(Yt_100,Ytmin1_100);
    rho_hat(i) = beta;
end
    
mean = mean(rho_hat)

I compared the values of the mean from these two codes and they are very different. Since I used random number generator, i get that it may be different but I do not expect the difference to be this large. For example, I get the mean from the first code to be 0.11 while the mean in the second code to be 0.45.

Thank you very much

Answer 1

[Cody LaFlamme]

In your first code block, when you call the sim_AR function, it seems you've switched alpha and rho accidentally. The provided function expects them to be passed in the opposite order. :)

[Walter Roberson]

Good catch!

[Ariza Gusti]

Thank you so much!

[Steven Lord]

FYI Ariza Gusti, if you're trying to debug code that calls the random number generator consider calling rng (with one or two input arguments) to set the random number generator to a known, fixed state before running your code. That eliminates differences in the sequence of random numbers generated by the rand* calls in your code as a potential cause of whatever you're investigating.