Group Abstract

Message Boards

WOLFRAM COMMUNITY

11.7K Views

12 Replies

6 Total Likes

View groups...

Follow this post

Share this post:

GROUPS:

Compiled function uses too much memory

Ulrich Utiger

Posted 6 years ago

This algorithm calculates the probability that S of n biological mtMRCA lineages survive after v generations. In the beginning (generation v=1), there are n mothers that have each randomly between 0 and d daughters. I have a fast analytical solution. But I need to check it with this algorithm. The uncompiled function is slow of course, but it does not need much memory. Unfortunately, memory usage of the compiled version is exponential for d > 2, n > 4 and V > 7 for a reason I don't understand. Can anyone help? ClearAll["Global`*"] Clear[distribution]; distribution = Compile[{{d, _Integer}, {n, _Integer}, {S, _Integer}, {z, _Integer}, \ {V, _Integer}}, Module[{w, r, v, s, c, K}, w = Table[0, {V}]; r = Table[0, {n}]; Do[ v = 1; s = 0; Do[ c = RandomInteger[{0, d}]; If[c != 0, s++; r[[k]] = c], {k, 1, n}]; If[s == S, w[[v]]++]; While[s != 0, v++; If[v == V, Break[]]; K = s; s = 0; Do[ c = Total[RandomInteger[{0, d}, r[[k]]]]; If[c != 0, s++; r[[k]] = c], {k, 1, K}]; If[s == S, w[[v]]++]], {z}]; Return[Table[{v, w[[v]]/z // N}, {v, 1, V - 1}]]], CompilationTarget -> "C"]

This algorithm calculates the probability that S of n biological mtMRCA lineages survive after v generations. In the beginning (generation v=1), there are n mothers that have each randomly between 0 and d daughters. I have a fast analytical solution. But I need to check it with this algorithm. The uncompiled function is slow of course, but it does not need much memory. Unfortunately, memory usage of the compiled version is exponential for d > 2, n > 4 and V > 7 for a reason I don't understand. Can anyone help?

ClearAll["Global`*"]
Clear[distribution];
distribution = 
Compile[{{d, _Integer}, {n, _Integer}, {S, _Integer}, {z, _Integer}, \
{V, _Integer}},
  Module[{w, r, v, s, c, K},
   w = Table[0, {V}]; r = Table[0, {n}];
   Do[
    v = 1; s = 0;
    Do[
     c = RandomInteger[{0, d}];
     If[c != 0, s++; r[[k]] = c],
     {k, 1, n}];
    If[s == S, w[[v]]++];
    While[s != 0,
     v++; If[v == V, Break[]];
     K = s; s = 0;
     Do[
      c = Total[RandomInteger[{0, d}, r[[k]]]];
      If[c != 0, s++; r[[k]] = c],
      {k, 1, K}];
     If[s == S, w[[v]]++]],
    {z}]; Return[Table[{v, w[[v]]/z // N}, {v, 1, V - 1}]]],
  CompilationTarget -> "C"]

POSTED BY: Ulrich Utiger

12 Replies

Sort By:

Ulrich Utiger

Posted 1 year ago

I finished my article about the probability of monogenism. I didn't have time to finish it 4 years ago but now it's done. Suggestions are welcome. Edit: the link should work now...

POSTED BY: Ulrich Utiger

Victor Kryukov

Victor Kryukov, PracticalWolf.com

Posted 1 year ago

POSTED BY: Victor Kryukov

Neil Singer

Neil Singer, AC Kinetics, Inc.

Posted 6 years ago

Ulrich, Glad this helps. I look forward to the article (although it is a bit outside of my area of expertise). Also, I would change the title of this post to make it more useful to others. Regards, Neil

POSTED BY: Neil Singer

Ulrich Utiger

Posted 6 years ago

POSTED BY: Ulrich Utiger

Neil Singer

Neil Singer, AC Kinetics, Inc.

Posted 6 years ago

POSTED BY: Neil Singer

Ulrich Utiger

Posted 6 years ago

Yes, I have VisualStudio installed. And to my surprise your code compiled flawlessly... Normally, when I did this in the past, it always refused some code that I needed to render C-compatible. The next step is now to find a probability distribution for the amount of daughters because in reality this is not a uniform probability as rendered with RandomInteger. A mother, for instance, that has 0 daughters is much less probable than if she has 3, in any case this is true for hunter-gatherers. I tested it with RandomVariate and the binomial distribution, which yields quite a different end-result than with RandomInteger.

POSTED BY: Ulrich Utiger

Ulrich Utiger

Posted 6 years ago

POSTED BY: Ulrich Utiger

Martijn Froeling

Martijn Froeling, University Medical Center Utrecht

Posted 6 years ago

As long as you have a C compiler installed on your system this is not difficult. Just check if one is availible using Needs["CCompilerDriver`"] DefaultCCompiler[] CCompilers[Full] I quickly checked and compiling with `CompilationTarget -> "C"` of the compiled function in my previous post makes it faster. (d3 = distributionF1C[3, 5, 3, 10000, 10]) // RepeatedTiming (d3 = distributionF1C2[3, 5, 3, 10000, 10]) // RepeatedTiming Out[17]= {0.12, {{1., 0.2657}, {2., 0.3424}, {3., 0.3544}, {4., 0.3604}, {5., 0.3613}, {6., 0.3611}, {7., 0.3586}, {8., 0.3576}, {9., 0.3574}}} Out[18]= {0.083, {{1., 0.2731}, {2., 0.3305}, {3., 0.3393}, {4., 0.3391}, {5., 0.3397}, {6., 0.3402}, {7., 0.3395}, {8., 0.3399}, {9., 0.3393}}}

As long as you have a C compiler installed on your system this is not difficult.

Just check if one is availible using

Needs["CCompilerDriver`"]
DefaultCCompiler[]
CCompilers[Full]

I quickly checked and compiling with CompilationTarget -> "C" of the compiled function in my previous post makes it faster.

(d3 = distributionF1C[3, 5, 3, 10000, 10]) // RepeatedTiming
(d3 = distributionF1C2[3, 5, 3, 10000, 10]) // RepeatedTiming

Out[17]= {0.12, {{1., 0.2657}, {2., 0.3424}, {3., 0.3544}, {4., 
   0.3604}, {5., 0.3613}, {6., 0.3611}, {7., 0.3586}, {8., 
   0.3576}, {9., 0.3574}}}

Out[18]= {0.083, {{1., 0.2731}, {2., 0.3305}, {3., 0.3393}, {4., 
   0.3391}, {5., 0.3397}, {6., 0.3402}, {7., 0.3395}, {8., 
   0.3399}, {9., 0.3393}}}

POSTED BY: Martijn Froeling