How to compare run-time of a model when run in MATLAB vs. OCTAVE?

help
#1

Thank you for writing us about your problem with using :octave: GNU Octave. To help us help you, please answer the following questions, if applicable.

Problem description

I have developed an agent-based model (ABM) which is optimized to run efficiently in MATLAB (R2020b). I am checking whether the model takes a similar amount of time when run in OCTAVE version 6.2.0. The ABM’s compute codes were adjusted with a couple of syntactical changes before running it in OCTAVE. Two changes were as follows: (1) discretize (MATLAB) was replaced with lookup (OCTAVE); and (2) MATLAB’s rng(seed) was replaced with OCTAVE’s rand(“seed”).

In MATLAB, the run-time for 4 replications of the ABM was 8036 seconds (Elapsed time is 8035.594101 seconds).

The first replication of the ABM did not finish running in OCTAVE by the time all the 4 replications finished running in MATLAB, even though it was started 60 minutes (clock time) before the start of the first replication in MATLAB.

How should I check run-time of the ABM in OCTAVE to compare it with that taken by the ABM when run in MATLAB?

My system

  • OS: Windows 10 Enterprise with 64-bit OS and processor: Intel(R) Xeon(R) CPU E5-2687W v3 @ 3.10GHz 3.10 GHz
  • Octave version: Version 6.2.0
  • Installation method: Downloaded and installed octave-6.2.0-w64-installer.exe" from Download
Reducing function copy overhead
#2

first, you should try the test on a trivially small data set. something that runs for hours on the faster system is probably not ideal for that. minutes or even seconds should be your target. you should run against something that still moves through most of the codepaths in the program, but will simply not take as long. then run that in both systems. there are many ways matlab is faster than octave. in particular, it has an effective Just-In-Time compiler that does optimizations on the fly in m-code that octave simply can’t do. It is why older methods focusing on vectorization are less necessary in Matlab than they used to be, but are still crucial for performance in Octave. Then I would make use of the profiler in both programs to see where the bottlenecks are, what parts are taking the most processing time. this can potentially be done on subfunctions, making the task a bit easier. Once you identify the bottlenecks, see if there are things you could do to streamline.

1 Like
#3

Dear @nrjank,

Thank you for your response. It’s very helpful.

My model’s computer codes are optimized for running in Matlab’s parallel environment. Using Matlab’s profiler I have profiled the codes several times. I am confident that any further attempt of profiling will result in zero or non-significant improvement in time taken by the model to finish a single simulation. Because the model is an agent-based model and it is solved using the Gillespie’s event driven algorithm, I believe vectorizing provides limited advantage. However, I am happy to profile the codes in Octave. I am not familiar about how to do it in Octave. Can you please point me to a correct documentation on how to profile computer codes in Octave? Thank you again.

Braj

#4

the octave profiler is much less sophisticated than the Matlab’s but it can be useful. It will give number of function calls and runtime spent on those function calls. it’s CLI only and doesn’t give you specific code line numbers, but the functions can be found in the Octave documentation:
https://octave.org/doc/latest/Profiling.html

the simple explanation is:

profile clear ## to start fresh if profile was run before
profile on
<run code>
profile off
profshow ##to get a simple summary
profexplore ##to get an environment where you can step down into each function call)

recommended only for code that runs to completion, so some trivial task compared to the one you mentioned at the top. you should be able to ctrl-c out of something running too long, then turn off the profiler and see what the times look like, but haven’t tried that and verified it works.

#5

Hi @nrjank,

Thanks again for your suggestion! I did profile a small part of my Agent-based model codes in OCTAVE 6.2.0 on my laptop with OS Microsoft Windows 10 Enterprise SP0 - 64 bit. The same part I also profiled in Matlab 2020b on the laptop. I have the following outcomes:

OCTAVE

>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 155.265 seconds.
>> profshow(2)
   #                                                                   Function Attr     Time (s)   Time (%)        Calls

-------------------------------------------------------------------------------------------------------------------------

  25 SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9           123.644      79.76        75018

   1                                         SyntheticIndividuals_for_ABM_debug            31.003      20.00            1

MATLAB

image
image1862×224 29.8 KB

I am not sure you can see the png file above for time spent during profiling of the part in Matlab. The difference is huge. I am also attaching the m-file and an input csv file in case you would want to check it for yourself.
SyntheticIndividuals_for_ABM_debug.m (7.1 KB)
input_data.csv (2.7 KB)

I have no idea at the moment at which place the program is taking so much time when run in OCTAVE. I would very much appreciate your insight.

Thanks.

#6

the png isn’t showing up, but did you try using profexplore as well? it lets you dive into each function call and see at least a little bit better how many times things are being called and where the time is being spent.

#7

Octave/matlab is not the most appropriate language for this type of problems. I would recommend Julia:

#8

Thank you, @nrjank!

Here are profiling results from another try.

Profiling results from OCTAVE:

Elapsed time is 159.834 seconds.
>> profshow(3)
   #                                                                   Function Attr     Time (s)   Time (%)        Calls
-------------------------------------------------------------------------------------------------------------------------
  25 SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9           125.720      78.79        75018
   1                                         SyntheticIndividuals_for_ABM_debug            33.412      20.94            1
  15                                                                   binary +             0.303       0.19      1350384
>> profexplore

Top
  1) SyntheticIndividuals_for_ABM_debug: 1 calls, 159.566 total, 33.412 self
  2) profshow: 1 calls, 0.004 total, 0.001 self
  3) profexplore: 1 calls, 0.000 total, 0.000 self

profexplore> 1

Top
  SyntheticIndividuals_for_ABM_debug: 1 calls, 159.566 total, 33.412 self
    1) SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9: 75018 calls, 126.095 total, 125.720 self
    2) binary +: 150036 calls, 0.053 total, 0.053 self
    3) SyntheticIndividuals_for_ABM_debug>SampleIndividualAges_From_AgeGroups: 40 calls, 0.002 total, 0.001 self
    4) SyntheticIndividuals_for_ABM_debug>get_Population_proportion_by_AgeGroup: 1 calls, 0.002 total, 0.001 self
    5) SyntheticIndividuals_for_ABM_debug>Create_Agents_of_a_community_v9: 1 calls, 0.001 total, 0.000 self
    6) toc: 1 calls, 0.000 total, 0.000 self
    7) prefix -: 40 calls, 0.000 total, 0.000 self
    8) sum: 1 calls, 0.000 total, 0.000 self
    9) ceil: 2 calls, 0.000 total, 0.000 self
    10) tic: 1 calls, 0.000 total, 0.000 self
    11) length: 4 calls, 0.000 total, 0.000 self
    12) binary /: 2 calls, 0.000 total, 0.000 self
    13) binary *: 2 calls, 0.000 total, 0.000 self

The profiling results from MATLAB are shown in a png file. I attach the png file this time.

Matlab_profiling_result
Matlab_profiling_result1912×741 95.8 KB

As you can see in the screenshot, the program finished in less than 0.3 sec in MATLAB. In OCTAVE, it took more than 150 sec. Which OCTAVE syntaxes are so time-expensive - is it rand?

#9

Thank you, @dasergatskov, for your suggestion! How easy it would be convert m-files into Julia? I am not familiar with Julia. (I have more than 6000 lines of computer codes of an agent-based model, developed and parameterized in MATLAB.

#10

in profexplore you can keep diving deeper. you need to dive into 1) again to look at

SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9: 75018 calls, 126.095 total, 125.720 self

and see what’s going on in there. so far this tells you nothing you don’t already know.

Octave can be very slow on looping structures and some other things, while Matlab has a JIT compiler that makes such issues mostly moot.

#11

Thank you, @nrjank, for your time and response.

I profiled the same part of my Agent-based model codes
SyntheticIndividuals_for_ABM_debug.m (7.1 KB)
in OCTAVE once again. The results are as follows:

>> profile on
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 170.842 seconds.
>> profshow(3)
   #                                                                   Function Attr     Time (s)   Time (%)        Calls
-------------------------------------------------------------------------------------------------------------------------
  25 SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9           136.491      80.03        75018
   1                                         SyntheticIndividuals_for_ABM_debug            33.608      19.70            1
  15                                                                   binary +             0.314       0.18      1350384

>> profexplore

Top
  1) SyntheticIndividuals_for_ABM_debug: 1 calls, 170.556 total, 33.608 self
  2) profshow: 1 calls, 0.004 total, 0.001 self
  3) profexplore: 1 calls, 0.000 total, 0.000 self

profexplore> 1

Top
  SyntheticIndividuals_for_ABM_debug: 1 calls, 170.556 total, 33.608 self
    1) SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9: 75018 calls, 136.881 total, 136.491 self
    2) binary +: 150036 calls, 0.057 total, 0.057 self
    3) SyntheticIndividuals_for_ABM_debug>get_Population_proportion_by_AgeGroup: 1 calls, 0.006 total, 0.001 self
    4) SyntheticIndividuals_for_ABM_debug>SampleIndividualAges_From_AgeGroups: 40 calls, 0.002 total, 0.001 self
    5) SyntheticIndividuals_for_ABM_debug>Create_Agents_of_a_community_v9: 1 calls, 0.002 total, 0.000 self
    6) toc: 1 calls, 0.000 total, 0.000 self
    7) prefix -: 40 calls, 0.000 total, 0.000 self
    8) sum: 1 calls, 0.000 total, 0.000 self
    9) ceil: 2 calls, 0.000 total, 0.000 self
    10) tic: 1 calls, 0.000 total, 0.000 self
    11) binary *: 2 calls, 0.000 total, 0.000 self
    12) length: 4 calls, 0.000 total, 0.000 self
    13) binary /: 2 calls, 0.000 total, 0.000 self

profexplore> 1

Top
  SyntheticIndividuals_for_ABM_debug: 1 calls, 170.556 total, 33.608 self
    SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9: 75018 calls, 136.881 total, 136.491 self
      1) binary +: 1200288 calls, 0.256 total, 0.256 self
      2) NaN: 75018 calls, 0.134 total, 0.134 self

profexplore> 1

Top
  SyntheticIndividuals_for_ABM_debug: 1 calls, 170.556 total, 33.608 self
    SyntheticIndividuals_for_ABM_debug>Add_an_individual_to_CommunityAgents_v9: 75018 calls, 136.881 total, 136.491 self
      binary +: 1200288 calls, 0.256 total, 0.256 self

As you can see, I do not seem to get any more deeper with profexplore. So, now it is clear to me that filling out (by calling Add_an_individual_to_CommunityAgents_v9) a total of 75018 rows in the data matrix AgentsCOMM (see the code, attached again) is more (>650x) time consuming than what (<0.2 sec) this process takes to complete in MATLAB. So, you are right. The code’s execution may not be sped up in OCTAVE while preserving the code’s structure. Is this a fair assessment?

#12

ok, i finally dug into the code.

at a first glance it looks like that function adds a record to the existing array. it appears all entries are numeric, and the item simply takes a [n x 15] array and outputs a [n+1 x 15 array] with the new info added.

If I’m reading this correctly, I can see why this would be slow. you are dynamically resizing an array. it has to allocate the new array and copy the old values ta each iteration. Matlab has made specific improvements it the past to minimize these penalties, and i suspect that matlab’s JIT compiler may be helping as well, but this used to be a big issue in slowing down Matlab code.

the best advise I can give is look at AgentsCOMM and see if there’s any way to preallocate the array size once you know how many times you’ll call the expansion function.

Edit to add - seeing as how you’re only calling that function in one place, and it’s in a for loop where you already know how many times you’ll call it, this looks like a plausible fix. it should run faster in both programs as well. there are quite a few octave and matlab resources out there on vectorizing for loops for improved performance. I suspect your program could take such and approach to eliminate the top level loops.

#13

This is an easy fix, and easy to test as well. It should make a difference.

Dynamically growing an array is stressful not so much on Octave, but on the underlying operating system libraries which have to manage the heap. Matlab also relies on the underlying OS for basic memory management so runtime there should also improve.

#14

Thank you, @nrjank, for taking time to look at my code. I really appreciate it. But I am not sure how my preallocation of the AgentsCOMM was missed.

In the code that I shared in my previous response, the preallocation is done in a function Create_Agents_of_a_community_v9(sum(PopulationSize)) by the following statements:

PopulationSize=75000;

AgentsCOMM=Create_Agents_of_a_community_v9(sum(PopulationSize));

However, I noticed that the code adds 18 additional rows to the preallocated AgentsCOMM. I fixed this issue by adding extra 50 rows at the time of preallocation as shown below.

PopulationSize=75000;    
AgentsCOMM=Create_Agents_of_a_community_v9(sum(PopulationSize)+50);

For ease, I provide the function Create_Agents_of_a_community_v9 here, which is defined as follows:

function AgentsCOMM = Create_Agents_of_a_community_v9( PopulationSize )
%% Preallocation and Inititalization
% preallocate the array
AgentCharacteristics={'ID';'DateOfBirth';'Gender';'LifeStatus';'Race';...
    'HospitalAdmissionTimes';'HospitalDischargeTimes';...
    'AbxHistory';'HealthStatus';'CurrentlyInHospital_or_Not';...
    'Which_HCF';'Which_ward';'CommunityID';'ClinicalTest';...
    'TestResultTime'};
AgentsCOMM=zeros(ceil(PopulationSize),length(AgentCharacteristics));
end

Even after the fix, elapsed time in completing the run in OCTAVE (without profile on) did not improve any significantly. Here are four examples.

>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 169.144 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 155.226 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 153.169 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 177.961 seconds.

Elapsed time for four equivalent runs in MATLAB are as follows:

>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 0.038057 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 0.049095 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 0.023833 seconds.
>> SyntheticIndividuals_for_ABM_debug();
Elapsed time is 0.025742 seconds.

These can be easily verified by anyone who has OCTAVE version 6.2.0 and/or MATLAB R2020b installed on their computers.

#15

I do not have Matlab and cannot do a side by side comparison. I can see that a lot of the time in Octave seems to be spent copying the AgentsCOMM matrix back and forth multiple times. I have a hypothesis as to why. Octave often does not make a fresh copy of a variable that is passed as an argument to a function unless that variable is modified inside the function. In your function to add an agent, you are using the same variable AgentsCOMM as input and output, likely causing extra copies to be made of the whole cell array. Try two things: first rename the variables in the function so that the input is not getting changed with the addition of a row. That may or may not speed it up. Second, eliminate that function to add a new row and absorb its function body inside the for loop where it is called. Namely, add row i+1 in the for loop directly not through a function call. You are already preallocating.

I know from past experience that Matlab is good at compiling functions and keeping them as compiled objects, so it has most likely eliminated the extra copying and is updating the array in place.

Edit: I have a suspicion that it’s scaling quadratically because of the array copying in Octave. Try running the Octave version with 5K, 10K, 20K etc agents. Do you expect linear time or quadratic time?

#16

As I expected, the change I suggested to eliminate that add function and move its code inside the location where it is called has sped it up by a lot. I deleted the function Add_an_individual_to_CommunityAgents_v9 completely and changed the main function’s for-loop as follows:

            iii=0;
            for ii=1:SummedAgeGrPopul_gender(iAgeGr,genderX)
                iii=iii+1;
                
                AgentsCOMM(ID+1,1)=ID+1;     % Identification number
                AgentsCOMM(ID+1,2)=SampledAges(iii);      % Death of birth
                AgentsCOMM(ID+1,3)=genderX;   % 1: Male or 2: Female
                AgentsCOMM(ID+1,4)=1;        % 1: alive 0: dead
                AgentsCOMM(ID+1,5)=nan;      % At the moment Race is NaN
                AgentsCOMM(ID+1,6)=0;        % Time of hospitalization as patient
                AgentsCOMM(ID+1,7)=0;        % Time of discharge from hospital as patient
                AgentsCOMM(ID+1,8)=0;        % Antibiotix history (0, 1, 2, or 3);
                AgentsCOMM(ID+1,9)=1;        % initialized with the susceptible status 
                AgentsCOMM(ID+1,10)=0;       % 0 currently NOT in a hospital; 1 in hospital 
                AgentsCOMM(ID+1,11)=0;       % ID of HCF person most recently admitted
                AgentsCOMM(ID+1,12)=0;       % ID of ward of the HCF
                AgentsCOMM(ID+1,13)=iComm;  % ID of the community a person lives in
                AgentsCOMM(ID+1,14)=0;     % Clinical Test result (the most recent)
                AgentsCOMM(ID+1,15)=0;     % Time of the most recent test result

                ID=ID+1;
            end

I renamed that version SyntheticIndividuals_for_ABM_TEST.m. Now see the performance. Each of the following is run with a fresh Octave instance:

Original:

octave:1> SyntheticIndividuals_for_ABM_debug
Elapsed time is 34.7479 seconds.

Improved:

octave:1> SyntheticIndividuals_for_ABM_TEST
Elapsed time is 1.71566 seconds.

That function caused unnecessary copying and was getting in the way. You can see why it becomes quadratic instead of linear if written like that. The entire array AgentsCOMM had to get copied from the main function SyntheticIndividuals_for_ABM_debug into the helper function Add_an_individual_to_CommunityAgents_v9, and it got copied back. So every addition of a row caused the whole matrix to be copied back and forth, meaning that the addition of N rows was causing the copying of (2 * N^2) rows’-worth of data.

Take-home lesson: Avoid function calls to grow an array like that. Inline it instead.

Technique lesson: If you see a slow-down that bad, experiment with the size of the problem, like changing 75K to 10K, 20K, etc. If you can see quadratic behavior for what should be mostly linear, you can expect that an entire array is being copied back and forth for some reason, so it’ll help localize the problem.

Edit:

I get into the time analysis in more detail here.

This is the original code behavior:

octave:17> pos = 0; x = y = [];
octave:18> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 10K
Elapsed time is 0.550006 seconds.
octave:19> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 20K
Elapsed time is 1.69111 seconds.
octave:20> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 30K
Elapsed time is 3.27898 seconds.
octave:21> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 40K
Elapsed time is 5.56754 seconds.
octave:22> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 50K
Elapsed time is 9.40509 seconds.
octave:23> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 60K
Elapsed time is 16.3135 seconds.
octave:24> SyntheticIndividuals_for_ABM_debug ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 70K
Elapsed time is 26.7577 seconds.
octave:25> [x; y]
ans =

    1.0000    2.0000    3.0000    4.0000    5.0000    6.0000    7.0000
    0.5500    1.6911    3.2795    5.5683    9.4059   16.3146   26.7588

octave:26> p = polyfit (x, y, 2), z = polyval (p, x); [x; y; z]
p =

   0.9073  -3.1871   3.6833

ans =

    1.0000    2.0000    3.0000    4.0000    5.0000    6.0000    7.0000
    0.5500    1.6911    3.2795    5.5683    9.4059   16.3146   26.7588
    1.4035    0.9384    2.2878    5.4519   10.4306   17.2240   25.8320

octave:27> p = polyfit (x, y, 3), z = polyval (p, x); [x; y; z]
p =

   0.1516  -0.9124   3.0301  -1.7757

ans =

    1.0000    2.0000    3.0000    4.0000    5.0000    6.0000    7.0000
    0.5500    1.6911    3.2795    5.5683    9.4059   16.3146   26.7588
    0.4937    1.8482    3.1977    5.4519    9.5208   16.3142   26.7419

As you can see, the time performance for the original is strongly quadratic and weakly cubic, so this is an indication that something weird is going on.

Here is the same analysis for the improved version, just from the change of inlining the code of that add_individual_to_group function:

octave:30> pos = 0; x = y = [];
octave:31> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 10K
Elapsed time is 0.225497 seconds.
octave:32> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 20K
Elapsed time is 0.445472 seconds.
octave:33> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 30K
Elapsed time is 0.66836 seconds.
octave:34> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 40K
Elapsed time is 0.889586 seconds.
octave:35> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 50K
Elapsed time is 1.11439 seconds.
octave:36> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 60K
Elapsed time is 1.35935 seconds.
octave:37> SyntheticIndividuals_for_ABM_TEST ; pos += 1; x(pos) = pos; y(pos) = toc;   ## 70K
Elapsed time is 1.56989 seconds.
octave:38> [x; y]
ans =

   1.0000   2.0000   3.0000   4.0000   5.0000   6.0000   7.0000
   0.2255   0.4455   0.6684   0.8896   1.1152   1.3603   1.5709

octave:39> p = polyfit (x, y, 2), z = polyval (p, x); [x; y; z]
p =

   8.6829e-04   2.1851e-01   5.1089e-03

ans =

   1.0000   2.0000   3.0000   4.0000   5.0000   6.0000   7.0000
   0.2255   0.4455   0.6684   0.8896   1.1152   1.3603   1.5709
   0.2245   0.4456   0.6684   0.8930   1.1194   1.3474   1.5772

octave:40> p = polyfit (x, y, 1), z = polyval (p, x); [x; y; z]
p =

   2.2545e-01  -5.3106e-03

ans =

   1.0000   2.0000   3.0000   4.0000   5.0000   6.0000   7.0000
   0.2255   0.4455   0.6684   0.8896   1.1152   1.3603   1.5709
   0.2201   0.4456   0.6710   0.8965   1.1220   1.3474   1.5729

It is strongly linear and the quadratic coefficient is tiny, which is much more plausible behavior.

All of the above was done on a single thread of a single core, no parallel computation or funny business.

1 Like
#17

Forgive me for coming back to this topic, but it has been bothering me with how it is written.

@SinghBK : If you eliminate the for-loop that has iii in it entirely and vectorize that loop as follows:

34,40c33,52
<             iii=0;
<             for ii=1:SummedAgeGrPopul_gender(iAgeGr,genderX)
<                 iii=iii+1;
<                 AgentsCOMM=Add_an_individual_to_CommunityAgents_v9(ID,...
<                     AgentsCOMM,genderX,SampledAges(iii),iComm);
<                 ID=ID+1;
<             end
---
> 
>             ii = (1:SummedAgeGrPopul_gender(iAgeGr,genderX))';
>             
>             AgentsCOMM(ID+ii,1) = ID + ii;  % Identification number
>             AgentsCOMM(ID+ii,2)=SampledAges(ii);      % Death of birth
>             AgentsCOMM(ID+ii,3)=genderX;   % 1: Male or 2: Female
>             AgentsCOMM(ID+ii,4)=1;        % 1: alive 0: dead
>             AgentsCOMM(ID+ii,5)=nan;      % At the moment Race is NaN
>             AgentsCOMM(ID+ii,6)=0;        % Time of hospitalization as patient
>             AgentsCOMM(ID+ii,7)=0;        % Time of discharge from hospital as patient
>             AgentsCOMM(ID+ii,8)=0;        % Antibiotix history (0, 1, 2, or 3);
>             AgentsCOMM(ID+ii,9)=1;        % initialized with the susceptible status 
>             AgentsCOMM(ID+ii,10)=0;       % 0 currently NOT in a hospital; 1 in hospital 
>             AgentsCOMM(ID+ii,11)=0;       % ID of HCF person most recently admitted
>             AgentsCOMM(ID+ii,12)=0;       % ID of ward of the HCF
>             AgentsCOMM(ID+ii,13)=iComm;  % ID of the community a person lives in
>             AgentsCOMM(ID+ii,14)=0;     % Clinical Test result (the most recent)
>             AgentsCOMM(ID+ii,15)=0;     % Time of the most recent test result
>             
>             ID += numel(ii);

Then not only do you get the same AgentsCOMM as before, but the timing performace is comparable to Matlab. These are the timings for a population of 75K:

Original code as submitted by @SinghBK :

octave:1> SyntheticIndividuals_for_ABM_debug
Elapsed time is 35.1406 seconds.

Yesterday’s improvement from me to eliminate a function and a lot of copy-by-value:

octave:1> SyntheticIndividuals_for_ABM_TEST
Elapsed time is 1.75943 seconds.

Today’s optimization from me to vectorize and remove that loop completely:

octave:1> SyntheticIndividuals_for_ABM_Opt
Elapsed time is 0.011699 seconds.

That’s a full 3000X faster.

Did you have some really compelling reason for writing it with multiple for-loops and function calls? It only seemed to cripple the performance if written that way. Does Matlab encourage you to write loops or did it invisibly optimize it away so the user never knew there was a problem?

Please check the results and adapt to your needs, but do avoid unnecessary loops and function calls when you can vectorize or broadcast arrays.