Saving variables from functions to workspace and making the code recognize the function

I'm trying to implement a logistic growth stochastic simulation, with my code being split into 4 functions as shown below. I'm getting one error in the code saying that function simulation might be unused. I've tried calling the function in the usual way as:

[n_vec]=simulation();

But that gives an 'undefined function or variable 'simulation' error. Any tips on how to successfully run the function and save the resulting n_vec to the workspace? I'd also like to save variables like br or dr, but I assume I can apply the same method as would be suggested for n_vec to do it on my own.

The code is:

%n-population size
%b- per capita birth rate
%d- per capita death rate at low density
%m- controls strength of density dependence 
%differential equation for a stochastic analogue simulation:
%dn/dt= (b-d)*n +m*(n^2)
%carrying capacity K= (b-d)/m
%Stochastic population will grow and then fluctuate around a 
%value a bit lower than K, but at some time will go extinct
b=2;    d=1;    m=0.01;
params=[b d m];
function [br]=BirthRate(n,b)
    br=b*n;
end

function [dr]= DeathRate(n,d,m)
    dr=n*d+m*(n^2);
end

function [n]= timestep(n,delta_t)%params is specified in the function code
    %delta_t is the time step
    %params is a list containing population dynamic parameters,
    %split into three columns:1- b, per capita birth rate,
    %2- d, per capita death rate, 3- m strength of density
    %params[:,1] is per capita birth rate
    %params[:,2] is per capita death rate etc.
    %Births
    %mean number of births is br*delta_t, relised number is Poisson 
    %distributed with this mean
    disp(n)
    lambda= BirthRate(n,b)*delta_t;
    pd = poissrnd(lambda);
    n =n + pd;
    disp(n) 
    % deaths
    %for a survival process with probability per unit time of 
    %dying H, the probability of surviving delta_t is:
    %ps- probability of surviving=exp(-H*delta_t)
    %the number surviving is a binomial random number,
    %with number of trials= population size and probability that
    %each trial succeeds being the survival probability ps
    %for each individual, probability of death per unit time is:
    %pod=dr/n
    p_survive=exp((-DeathRate(n, d, m)*delta_t)/n); %check if brackets are corrects
    disp(p_survive)

      n=binornd(n,p_survive); %binomial random numbers

end
%code to run the simulation 
function [n_vec] = simulation() %n and params are specified in function, dont have to be added
  n=10; %initial population size
  n_vec=zeros(1,1000);
  for i = 1:1:1000 
      n=timestep(n,0.02, params);
      n_vec(i)=n;%adds the new(current) population size to the vector data
  end
end

Answer 1

[Steven Lord]

From the documentation:

"MATLAB® program files can contain code for more than one function. In a function file, the first function in the file is called the main function. This function is visible to functions in other files, or you can call it from the command line. Additional functions within the file are called local functions, and they can occur in any order after the main function. Local functions are only visible to other functions in the same file."

Your simulation function is a local function in this script file, and as such is not visible to functions outside this file. You have a couple options if you want to make it available outside, some of which are:

1. Make it its own main function in its own function file.
2. Have your script create a function handle for each local function and call the local functions using the handles. The second example on the documentation page for the localfunctions function shows an example of this.
3. Convert your script into a class file with simulation as a method (either a Static method that you call using the name of your class or an instance method that you must call with an instance of the class as one of the inputs.)

The first two approaches are easier. The last would require more work.

[Charles]

Thanks, I decided to go with option 2 and successfully created local handles. However now that I have access to them, I'm having some issues with actually running them. No matter what I try, I keep getting errors about insufficient inputs or undefined variables (despite them being constants in the main function). My current code is as below:

%n-population size
%b- per capita birth rate
%d- per capita death rate at low density
%m- controls strength of density dependence 
%differential equation for a stochastic analogue simulation:
%dn/dt= (b-d)*n +m*(n^2)
%carrying capacity K= (b-d)/m
%Stochastic population will grow and then fluctuate around a 
%value a bit lower than K, but at some time will go extinct
b=2;    d=1;    m=0.01;
params=[b d m];
delta_t=0.2;
fh=localfunctions; %store all local functions
n_vector = fh{4}(n,delta_t,params);
Carrying_Capacity=(b-d)/m;
function [br]=BirthRate(n,b)
    br=b*n;
end

function [dr]= DeathRate(n,d,m)
    dr=n*d+m*(n^2);
end

function [n]= timestep(n,delta_t,params)%params is specified in the function code
    %delta_t is the time step
    %params is a list containing population dynamic parameters,
    %mean number of births is br*delta_t, relised number is Poisson 
    %distributed with this mean
    disp(n)
    lambda= BirthRate(n,b)*delta_t;
    pd = poissrnd(lambda);
    n =n + pd;
    disp(n) 
    % deaths
    %for a survival process with probability per unit time of 
    %dying H, the probability of surviving delta_t is:
    %ps- probability of surviving=exp(-H*delta_t)
    %pod=dr/n
    p_survive=exp((-DeathRate(n, d, m)*delta_t)/n); %check if brackets are corrects
    disp(p_survive)
    n=binornd(n,p_survive); %binomial random numbers
end
%code to run the simulation, maybe not necessary as a seperate function?    
function [n_vec] = simulation(n,delta_t,params) %n and params are specified in function, dont have to be added
  %n=10; %initial population size
  n_vec=zeros(1,1000);
  for i = 1:1:1000 
      n=timestep(n,delta_t,params);
      n_vec(i)=n;%adds the new(current) population size to the vector data
  end
end

Any suggestions on how I could fix it? Looking at the example in the localfunctions page, when the handle is called specific values are provided, however in my code I can't provide n as a constant, as it's changing with time according to the timestep function.

[Steven Lord]

Where did you define the variable n that you pass as the first input to the function handle?

Can you show the full text of the error message? Often the exact wording provides details that help narrow down the problem, so please copy and paste rather than paraphrasing.

[Charles]

The variable n is defined first in the timestep function and is equal to it. My goal in the n_vec equation is to have an array of values that are: n value added to n_vec's i'th element (for the loop i = 1:1:1000).

And yes, of course. With the code identical to the one provided above, i get the error:

Undefined function or variable 'n'.

Error in stochastic_simulation_finaldraft (line 14) n_vector = fh{4}(n,delta_t,params);

I can also provide a reference on how the n_vec array should look like if needed. I have a working code developed in R and my goal is to successfully translate that into MATLAB code.

[Steven Lord]

You will need to call the timestep function to compute the value of n then pass that value into the simulation function. Calling simulation doesn't automatically call all the code prior to it in the file; if you want the timestep function to execute, you need to call it (or call a function that calls it, or call a function that calls a function that calls it, ... etc.)

n = fh{3}(... % Add the input arguments needed to call timestep
n_vector = fh{4}(n, delta_t, params);

[Charles]

Thank you, that did the trick once I adjusted the list of inputs and manually put constants in each of the functions, instead of just the main file.