Survival Probability for Random Walks












5














The Survival Probability for a walker starting at the origin is defined as the probability that the walker stays positive through n steps. Thanks to the Sparre-Andersen Theorem I know this PDF is given by



Plot[Binomial[2 n, n]*2^(-2 n), {n, 0, 100}]


However, I want to validate this empirically.



My attempt to validate this for n=100:



FoldList[
If[#2 < 0, 0, #1 + #2] &,
Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


I wantFoldList to stop if #2 < 0 evaluates to True, not just substitute in 0.










share|improve this question
























  • Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
    – MikeY
    2 days ago
















5














The Survival Probability for a walker starting at the origin is defined as the probability that the walker stays positive through n steps. Thanks to the Sparre-Andersen Theorem I know this PDF is given by



Plot[Binomial[2 n, n]*2^(-2 n), {n, 0, 100}]


However, I want to validate this empirically.



My attempt to validate this for n=100:



FoldList[
If[#2 < 0, 0, #1 + #2] &,
Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


I wantFoldList to stop if #2 < 0 evaluates to True, not just substitute in 0.










share|improve this question
























  • Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
    – MikeY
    2 days ago














5












5








5


2





The Survival Probability for a walker starting at the origin is defined as the probability that the walker stays positive through n steps. Thanks to the Sparre-Andersen Theorem I know this PDF is given by



Plot[Binomial[2 n, n]*2^(-2 n), {n, 0, 100}]


However, I want to validate this empirically.



My attempt to validate this for n=100:



FoldList[
If[#2 < 0, 0, #1 + #2] &,
Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


I wantFoldList to stop if #2 < 0 evaluates to True, not just substitute in 0.










share|improve this question















The Survival Probability for a walker starting at the origin is defined as the probability that the walker stays positive through n steps. Thanks to the Sparre-Andersen Theorem I know this PDF is given by



Plot[Binomial[2 n, n]*2^(-2 n), {n, 0, 100}]


However, I want to validate this empirically.



My attempt to validate this for n=100:



FoldList[
If[#2 < 0, 0, #1 + #2] &,
Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


I wantFoldList to stop if #2 < 0 evaluates to True, not just substitute in 0.







functions probability-or-statistics random distributions random-process






share|improve this question















share|improve this question













share|improve this question




share|improve this question








edited yesterday









yosimitsu kodanuri

437312




437312










asked 2 days ago









WillWill

1004




1004












  • Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
    – MikeY
    2 days ago


















  • Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
    – MikeY
    2 days ago
















Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
– MikeY
2 days ago




Will, are you attempting to empirically show that the probability for survival when n=100 is Binomial[2 (100), (100)]*2^(-2 (100))? So repeatedly run, and count the times you survive through 100 steps? If so, are you trying to "While" out of the FoldList to save CPU cycles? Not clear to me...
– MikeY
2 days ago










4 Answers
4






active

oldest

votes


















8














Something seems odd to me about your code. You are summing twice, once with Accumulate and once with FoldList. If this is really what you want then you could use:



SeedRandom[26]
sum = Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0];

TakeWhile[sum, NonNegative] // Accumulate



8

{0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964}



This is equivalent to your FoldList construct up to the appropriate point:



FoldList[If[#2 < 0, 0, #1 + #2] &, sum]



{0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964, 0, ...



Perhaps you meant to only sum once. In that case TakeWhile[sum, NonNegative] is a direct solution but also sub-optimal as it does not provide early exit behavior, which I suspect is what you're actually after here. It is not clear to me if you need the cumulative sum (walk) itself or only its length; if the latter consider this:



SeedRandom[26]
dist = RandomVariate[NormalDistribution[0, 1], 100];

Module[{i = 0},
Fold[If[# < 0, Return[i, Fold], i++; # + #2] &, 0, dist]
]



8






share|improve this answer





















  • The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
    – Carl Lange
    yesterday



















5














We can do this using an implementation of FoldWhileList.



First, implement FoldWhileList using this great answer.



FoldWhileList[f_, test_, start_, secargs_List] := 
Module[{tag},
If[# === {}, {start}, Prepend[First@#, start]] &@
Reap[Fold[If[test[##], Sow[f[##], tag], Return[Null, Fold]] &,
start, secargs], _, #2 &][[2]]]


Now we simply run this using the test #2 >= 0 (note that the implementation of NestWhile breaks when test stops evaluating True - our implementation of FoldWhileList also does this, therefore we invert the test you originally used.



FoldWhileList[Plus, #2 >= 0 &, 0, 
Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


We can now estimate your PDF:



pdf estimate



and overlay it over the original plot also:



overlaid plots



which doesn't seem like a great match - perhaps there's an issue with your original code, as this answer surmises.






share|improve this answer































    4














    It seems to me that this is a problem to which Catch and Throw can be usefully applied.



    SeedRandom[1];
    Module[{result = {0}, s},
    Catch[
    Fold[
    If[#2 < 0, Throw[Null], result = {result, s = #1 + #2}; s] &,
    0,
    Accumulate[RandomVariate[NormalDistribution[0, 1], 100]]]];
    result // Flatten]


    result






    share|improve this answer































      2














      How about the following brute force approach:



      n = 100;
      SeedRandom[12345];
      nsim = 1000000;
      Total[Table[If[Min[Accumulate[RandomVariate[NormalDistribution[0, 1], n]]] <= 0, 0, 1],
      {i, nsim}]]/nsim // N
      (* 0.056092 *)
      Binomial[2 n, n] 2^(-2 n) // N
      (* 0.0563485 *)


      To get all of the values from 1 to 100 "simultaneously"...



      SeedRandom[12345];
      nsim = 100000;
      n = 100;
      z = ConstantArray[0, n];
      Do[
      x = Accumulate[RandomVariate[NormalDistribution[0, 1], n]];
      i = Flatten[Position[x, _?NonPositive]];
      If[Length[i] > 0, If[i[[1]] > 1, z[[1 ;; i[[1]] - 1]] = z[[1 ;; i[[1]] - 1]] + 1], z = z + 1],
      {j, nsim}]
      z = z/nsim;
      ListPlot[{z, Table[Binomial[2 j, j] 2^(-2 j), {j, n}]}, PlotRange -> All, ImageSize -> Large]


      Simulation and exact formula






      share|improve this answer























        Your Answer





        StackExchange.ifUsing("editor", function () {
        return StackExchange.using("mathjaxEditing", function () {
        StackExchange.MarkdownEditor.creationCallbacks.add(function (editor, postfix) {
        StackExchange.mathjaxEditing.prepareWmdForMathJax(editor, postfix, [["$", "$"], ["\\(","\\)"]]);
        });
        });
        }, "mathjax-editing");

        StackExchange.ready(function() {
        var channelOptions = {
        tags: "".split(" "),
        id: "387"
        };
        initTagRenderer("".split(" "), "".split(" "), channelOptions);

        StackExchange.using("externalEditor", function() {
        // Have to fire editor after snippets, if snippets enabled
        if (StackExchange.settings.snippets.snippetsEnabled) {
        StackExchange.using("snippets", function() {
        createEditor();
        });
        }
        else {
        createEditor();
        }
        });

        function createEditor() {
        StackExchange.prepareEditor({
        heartbeatType: 'answer',
        autoActivateHeartbeat: false,
        convertImagesToLinks: false,
        noModals: true,
        showLowRepImageUploadWarning: true,
        reputationToPostImages: null,
        bindNavPrevention: true,
        postfix: "",
        imageUploader: {
        brandingHtml: "Powered by u003ca class="icon-imgur-white" href="https://imgur.com/"u003eu003c/au003e",
        contentPolicyHtml: "User contributions licensed under u003ca href="https://creativecommons.org/licenses/by-sa/3.0/"u003ecc by-sa 3.0 with attribution requiredu003c/au003e u003ca href="https://stackoverflow.com/legal/content-policy"u003e(content policy)u003c/au003e",
        allowUrls: true
        },
        onDemand: true,
        discardSelector: ".discard-answer"
        ,immediatelyShowMarkdownHelp:true
        });


        }
        });














        draft saved

        draft discarded


















        StackExchange.ready(
        function () {
        StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmathematica.stackexchange.com%2fquestions%2f189069%2fsurvival-probability-for-random-walks%23new-answer', 'question_page');
        }
        );

        Post as a guest















        Required, but never shown

























        4 Answers
        4






        active

        oldest

        votes








        4 Answers
        4






        active

        oldest

        votes









        active

        oldest

        votes






        active

        oldest

        votes









        8














        Something seems odd to me about your code. You are summing twice, once with Accumulate and once with FoldList. If this is really what you want then you could use:



        SeedRandom[26]
        sum = Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0];

        TakeWhile[sum, NonNegative] // Accumulate



        8

        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964}



        This is equivalent to your FoldList construct up to the appropriate point:



        FoldList[If[#2 < 0, 0, #1 + #2] &, sum]



        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964, 0, ...



        Perhaps you meant to only sum once. In that case TakeWhile[sum, NonNegative] is a direct solution but also sub-optimal as it does not provide early exit behavior, which I suspect is what you're actually after here. It is not clear to me if you need the cumulative sum (walk) itself or only its length; if the latter consider this:



        SeedRandom[26]
        dist = RandomVariate[NormalDistribution[0, 1], 100];

        Module[{i = 0},
        Fold[If[# < 0, Return[i, Fold], i++; # + #2] &, 0, dist]
        ]



        8






        share|improve this answer





















        • The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
          – Carl Lange
          yesterday
















        8














        Something seems odd to me about your code. You are summing twice, once with Accumulate and once with FoldList. If this is really what you want then you could use:



        SeedRandom[26]
        sum = Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0];

        TakeWhile[sum, NonNegative] // Accumulate



        8

        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964}



        This is equivalent to your FoldList construct up to the appropriate point:



        FoldList[If[#2 < 0, 0, #1 + #2] &, sum]



        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964, 0, ...



        Perhaps you meant to only sum once. In that case TakeWhile[sum, NonNegative] is a direct solution but also sub-optimal as it does not provide early exit behavior, which I suspect is what you're actually after here. It is not clear to me if you need the cumulative sum (walk) itself or only its length; if the latter consider this:



        SeedRandom[26]
        dist = RandomVariate[NormalDistribution[0, 1], 100];

        Module[{i = 0},
        Fold[If[# < 0, Return[i, Fold], i++; # + #2] &, 0, dist]
        ]



        8






        share|improve this answer





















        • The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
          – Carl Lange
          yesterday














        8












        8








        8






        Something seems odd to me about your code. You are summing twice, once with Accumulate and once with FoldList. If this is really what you want then you could use:



        SeedRandom[26]
        sum = Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0];

        TakeWhile[sum, NonNegative] // Accumulate



        8

        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964}



        This is equivalent to your FoldList construct up to the appropriate point:



        FoldList[If[#2 < 0, 0, #1 + #2] &, sum]



        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964, 0, ...



        Perhaps you meant to only sum once. In that case TakeWhile[sum, NonNegative] is a direct solution but also sub-optimal as it does not provide early exit behavior, which I suspect is what you're actually after here. It is not clear to me if you need the cumulative sum (walk) itself or only its length; if the latter consider this:



        SeedRandom[26]
        dist = RandomVariate[NormalDistribution[0, 1], 100];

        Module[{i = 0},
        Fold[If[# < 0, Return[i, Fold], i++; # + #2] &, 0, dist]
        ]



        8






        share|improve this answer












        Something seems odd to me about your code. You are summing twice, once with Accumulate and once with FoldList. If this is really what you want then you could use:



        SeedRandom[26]
        sum = Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0];

        TakeWhile[sum, NonNegative] // Accumulate



        8

        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964}



        This is equivalent to your FoldList construct up to the appropriate point:



        FoldList[If[#2 < 0, 0, #1 + #2] &, sum]



        {0, 1.10708, 1.23211, 2.28173, 3.30295, 4.05759, 5.26123, 6.62964, 0, ...



        Perhaps you meant to only sum once. In that case TakeWhile[sum, NonNegative] is a direct solution but also sub-optimal as it does not provide early exit behavior, which I suspect is what you're actually after here. It is not clear to me if you need the cumulative sum (walk) itself or only its length; if the latter consider this:



        SeedRandom[26]
        dist = RandomVariate[NormalDistribution[0, 1], 100];

        Module[{i = 0},
        Fold[If[# < 0, Return[i, Fold], i++; # + #2] &, 0, dist]
        ]



        8







        share|improve this answer












        share|improve this answer



        share|improve this answer










        answered 2 days ago









        Mr.WizardMr.Wizard

        230k294741038




        230k294741038












        • The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
          – Carl Lange
          yesterday


















        • The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
          – Carl Lange
          yesterday
















        The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
        – Carl Lange
        yesterday




        The ListLinePlot of your bottom answer is much closer to the binomial than the one in mine. I think you have it right :)
        – Carl Lange
        yesterday











        5














        We can do this using an implementation of FoldWhileList.



        First, implement FoldWhileList using this great answer.



        FoldWhileList[f_, test_, start_, secargs_List] := 
        Module[{tag},
        If[# === {}, {start}, Prepend[First@#, start]] &@
        Reap[Fold[If[test[##], Sow[f[##], tag], Return[Null, Fold]] &,
        start, secargs], _, #2 &][[2]]]


        Now we simply run this using the test #2 >= 0 (note that the implementation of NestWhile breaks when test stops evaluating True - our implementation of FoldWhileList also does this, therefore we invert the test you originally used.



        FoldWhileList[Plus, #2 >= 0 &, 0, 
        Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


        We can now estimate your PDF:



        pdf estimate



        and overlay it over the original plot also:



        overlaid plots



        which doesn't seem like a great match - perhaps there's an issue with your original code, as this answer surmises.






        share|improve this answer




























          5














          We can do this using an implementation of FoldWhileList.



          First, implement FoldWhileList using this great answer.



          FoldWhileList[f_, test_, start_, secargs_List] := 
          Module[{tag},
          If[# === {}, {start}, Prepend[First@#, start]] &@
          Reap[Fold[If[test[##], Sow[f[##], tag], Return[Null, Fold]] &,
          start, secargs], _, #2 &][[2]]]


          Now we simply run this using the test #2 >= 0 (note that the implementation of NestWhile breaks when test stops evaluating True - our implementation of FoldWhileList also does this, therefore we invert the test you originally used.



          FoldWhileList[Plus, #2 >= 0 &, 0, 
          Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


          We can now estimate your PDF:



          pdf estimate



          and overlay it over the original plot also:



          overlaid plots



          which doesn't seem like a great match - perhaps there's an issue with your original code, as this answer surmises.






          share|improve this answer


























            5












            5








            5






            We can do this using an implementation of FoldWhileList.



            First, implement FoldWhileList using this great answer.



            FoldWhileList[f_, test_, start_, secargs_List] := 
            Module[{tag},
            If[# === {}, {start}, Prepend[First@#, start]] &@
            Reap[Fold[If[test[##], Sow[f[##], tag], Return[Null, Fold]] &,
            start, secargs], _, #2 &][[2]]]


            Now we simply run this using the test #2 >= 0 (note that the implementation of NestWhile breaks when test stops evaluating True - our implementation of FoldWhileList also does this, therefore we invert the test you originally used.



            FoldWhileList[Plus, #2 >= 0 &, 0, 
            Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


            We can now estimate your PDF:



            pdf estimate



            and overlay it over the original plot also:



            overlaid plots



            which doesn't seem like a great match - perhaps there's an issue with your original code, as this answer surmises.






            share|improve this answer














            We can do this using an implementation of FoldWhileList.



            First, implement FoldWhileList using this great answer.



            FoldWhileList[f_, test_, start_, secargs_List] := 
            Module[{tag},
            If[# === {}, {start}, Prepend[First@#, start]] &@
            Reap[Fold[If[test[##], Sow[f[##], tag], Return[Null, Fold]] &,
            start, secargs], _, #2 &][[2]]]


            Now we simply run this using the test #2 >= 0 (note that the implementation of NestWhile breaks when test stops evaluating True - our implementation of FoldWhileList also does this, therefore we invert the test you originally used.



            FoldWhileList[Plus, #2 >= 0 &, 0, 
            Prepend[Accumulate[RandomVariate[NormalDistribution[0, 1], 100]], 0]]


            We can now estimate your PDF:



            pdf estimate



            and overlay it over the original plot also:



            overlaid plots



            which doesn't seem like a great match - perhaps there's an issue with your original code, as this answer surmises.







            share|improve this answer














            share|improve this answer



            share|improve this answer








            edited yesterday

























            answered 2 days ago









            Carl LangeCarl Lange

            2,2911625




            2,2911625























                4














                It seems to me that this is a problem to which Catch and Throw can be usefully applied.



                SeedRandom[1];
                Module[{result = {0}, s},
                Catch[
                Fold[
                If[#2 < 0, Throw[Null], result = {result, s = #1 + #2}; s] &,
                0,
                Accumulate[RandomVariate[NormalDistribution[0, 1], 100]]]];
                result // Flatten]


                result






                share|improve this answer




























                  4














                  It seems to me that this is a problem to which Catch and Throw can be usefully applied.



                  SeedRandom[1];
                  Module[{result = {0}, s},
                  Catch[
                  Fold[
                  If[#2 < 0, Throw[Null], result = {result, s = #1 + #2}; s] &,
                  0,
                  Accumulate[RandomVariate[NormalDistribution[0, 1], 100]]]];
                  result // Flatten]


                  result






                  share|improve this answer


























                    4












                    4








                    4






                    It seems to me that this is a problem to which Catch and Throw can be usefully applied.



                    SeedRandom[1];
                    Module[{result = {0}, s},
                    Catch[
                    Fold[
                    If[#2 < 0, Throw[Null], result = {result, s = #1 + #2}; s] &,
                    0,
                    Accumulate[RandomVariate[NormalDistribution[0, 1], 100]]]];
                    result // Flatten]


                    result






                    share|improve this answer














                    It seems to me that this is a problem to which Catch and Throw can be usefully applied.



                    SeedRandom[1];
                    Module[{result = {0}, s},
                    Catch[
                    Fold[
                    If[#2 < 0, Throw[Null], result = {result, s = #1 + #2}; s] &,
                    0,
                    Accumulate[RandomVariate[NormalDistribution[0, 1], 100]]]];
                    result // Flatten]


                    result







                    share|improve this answer














                    share|improve this answer



                    share|improve this answer








                    edited 2 days ago

























                    answered 2 days ago









                    m_goldbergm_goldberg

                    84.5k872196




                    84.5k872196























                        2














                        How about the following brute force approach:



                        n = 100;
                        SeedRandom[12345];
                        nsim = 1000000;
                        Total[Table[If[Min[Accumulate[RandomVariate[NormalDistribution[0, 1], n]]] <= 0, 0, 1],
                        {i, nsim}]]/nsim // N
                        (* 0.056092 *)
                        Binomial[2 n, n] 2^(-2 n) // N
                        (* 0.0563485 *)


                        To get all of the values from 1 to 100 "simultaneously"...



                        SeedRandom[12345];
                        nsim = 100000;
                        n = 100;
                        z = ConstantArray[0, n];
                        Do[
                        x = Accumulate[RandomVariate[NormalDistribution[0, 1], n]];
                        i = Flatten[Position[x, _?NonPositive]];
                        If[Length[i] > 0, If[i[[1]] > 1, z[[1 ;; i[[1]] - 1]] = z[[1 ;; i[[1]] - 1]] + 1], z = z + 1],
                        {j, nsim}]
                        z = z/nsim;
                        ListPlot[{z, Table[Binomial[2 j, j] 2^(-2 j), {j, n}]}, PlotRange -> All, ImageSize -> Large]


                        Simulation and exact formula






                        share|improve this answer




























                          2














                          How about the following brute force approach:



                          n = 100;
                          SeedRandom[12345];
                          nsim = 1000000;
                          Total[Table[If[Min[Accumulate[RandomVariate[NormalDistribution[0, 1], n]]] <= 0, 0, 1],
                          {i, nsim}]]/nsim // N
                          (* 0.056092 *)
                          Binomial[2 n, n] 2^(-2 n) // N
                          (* 0.0563485 *)


                          To get all of the values from 1 to 100 "simultaneously"...



                          SeedRandom[12345];
                          nsim = 100000;
                          n = 100;
                          z = ConstantArray[0, n];
                          Do[
                          x = Accumulate[RandomVariate[NormalDistribution[0, 1], n]];
                          i = Flatten[Position[x, _?NonPositive]];
                          If[Length[i] > 0, If[i[[1]] > 1, z[[1 ;; i[[1]] - 1]] = z[[1 ;; i[[1]] - 1]] + 1], z = z + 1],
                          {j, nsim}]
                          z = z/nsim;
                          ListPlot[{z, Table[Binomial[2 j, j] 2^(-2 j), {j, n}]}, PlotRange -> All, ImageSize -> Large]


                          Simulation and exact formula






                          share|improve this answer


























                            2












                            2








                            2






                            How about the following brute force approach:



                            n = 100;
                            SeedRandom[12345];
                            nsim = 1000000;
                            Total[Table[If[Min[Accumulate[RandomVariate[NormalDistribution[0, 1], n]]] <= 0, 0, 1],
                            {i, nsim}]]/nsim // N
                            (* 0.056092 *)
                            Binomial[2 n, n] 2^(-2 n) // N
                            (* 0.0563485 *)


                            To get all of the values from 1 to 100 "simultaneously"...



                            SeedRandom[12345];
                            nsim = 100000;
                            n = 100;
                            z = ConstantArray[0, n];
                            Do[
                            x = Accumulate[RandomVariate[NormalDistribution[0, 1], n]];
                            i = Flatten[Position[x, _?NonPositive]];
                            If[Length[i] > 0, If[i[[1]] > 1, z[[1 ;; i[[1]] - 1]] = z[[1 ;; i[[1]] - 1]] + 1], z = z + 1],
                            {j, nsim}]
                            z = z/nsim;
                            ListPlot[{z, Table[Binomial[2 j, j] 2^(-2 j), {j, n}]}, PlotRange -> All, ImageSize -> Large]


                            Simulation and exact formula






                            share|improve this answer














                            How about the following brute force approach:



                            n = 100;
                            SeedRandom[12345];
                            nsim = 1000000;
                            Total[Table[If[Min[Accumulate[RandomVariate[NormalDistribution[0, 1], n]]] <= 0, 0, 1],
                            {i, nsim}]]/nsim // N
                            (* 0.056092 *)
                            Binomial[2 n, n] 2^(-2 n) // N
                            (* 0.0563485 *)


                            To get all of the values from 1 to 100 "simultaneously"...



                            SeedRandom[12345];
                            nsim = 100000;
                            n = 100;
                            z = ConstantArray[0, n];
                            Do[
                            x = Accumulate[RandomVariate[NormalDistribution[0, 1], n]];
                            i = Flatten[Position[x, _?NonPositive]];
                            If[Length[i] > 0, If[i[[1]] > 1, z[[1 ;; i[[1]] - 1]] = z[[1 ;; i[[1]] - 1]] + 1], z = z + 1],
                            {j, nsim}]
                            z = z/nsim;
                            ListPlot[{z, Table[Binomial[2 j, j] 2^(-2 j), {j, n}]}, PlotRange -> All, ImageSize -> Large]


                            Simulation and exact formula







                            share|improve this answer














                            share|improve this answer



                            share|improve this answer








                            edited yesterday

























                            answered 2 days ago









                            JimBJimB

                            17.1k12663




                            17.1k12663






























                                draft saved

                                draft discarded




















































                                Thanks for contributing an answer to Mathematica Stack Exchange!


                                • Please be sure to answer the question. Provide details and share your research!

                                But avoid



                                • Asking for help, clarification, or responding to other answers.

                                • Making statements based on opinion; back them up with references or personal experience.


                                Use MathJax to format equations. MathJax reference.


                                To learn more, see our tips on writing great answers.




                                draft saved


                                draft discarded














                                StackExchange.ready(
                                function () {
                                StackExchange.openid.initPostLogin('.new-post-login', 'https%3a%2f%2fmathematica.stackexchange.com%2fquestions%2f189069%2fsurvival-probability-for-random-walks%23new-answer', 'question_page');
                                }
                                );

                                Post as a guest















                                Required, but never shown





















































                                Required, but never shown














                                Required, but never shown












                                Required, but never shown







                                Required, but never shown

































                                Required, but never shown














                                Required, but never shown












                                Required, but never shown







                                Required, but never shown







                                Popular posts from this blog

                                An IMO inspired problem

                                Management

                                Investment