📈 Simulating Poisson process (part 2)

[This article was first published on Iegor Rudnytskyi, and kindly contributed to R-bloggers]. (You can report issue about the content on this page here)
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

In previous post we discussed two common methods of Poisson process simulation. The reason why this trivial problem was of my interest is the fact that this is simplification of a larger scale problem of a classical ruin process. Let me remind that I focus on an extenssion of Cramér–Lundberg model with positive jumps, that is:

where:

  • $u$: is an initial capital;
  • $c$ is a premium rate;
  • $N_1(t)$ and $N_2(t)$ are Poisson processes of capital injections and claims with rates $\lambda_1$ and $\lambda_2$, respectively;
  • $X_i$ and $Y_j$ are i.i.d. random variables modeling sizes of capital injections and claims, respectively.

We simplify this model to a bare minimum, which reflects the behaviour of Poisson processes, that is we set $u = 0$, $c = 0$. Further, we assume deterministic unit jumps: $X_i = 1$ and $Y_j = 1$. As result we have the following model:

Nice thing about this model is that we know exact distribution of $X(t)$. It is called Skellam distribution, which is covered in skellam package. Therefore, we can compare Monte-Carlo simulated estimates to their exact equivaletns.

Method 1

For Gerber-Shui function we need to simulate a path until the ruin. It means that the time horizon is not known a priori. As consequence, the algorithm that first simulates the number of jumps for a given time is not applicable. Therefore, the only possibility to simulate a path until the ruin is to exploit the fact that interarrival times of jumps are exponentially distributed. In case of only negative jumps the algorithm is quite simple: in while loop we add jumps to a path until the process is ruined (or any other stopping conditions, e.g. maximum iterations acheived).

Things are slightly more complicated when the model includes positive jumps. If both positive and negative jumps’ times were known to us, then it would be possible to sort them, and add to a path in ascending order, as in an illustration below.

However, we do not know the arrival times of jumps, as they should be simulated. Also it would be naive to add a jump of one type at an iteration and then of another type in the next iteration, because underlying Possion process can have different rates (and, therefore, there is no garantee that one type jump follows the other type). The approach I propose here is very similar to a playground game “tag”. We generate arrival jump’s time for both types. Then, for one that occurred earlier (A), we need to catch up with the opposite type’s (B) time, that is generate more jumps of type (A) until the time of the later type (B) is achieved.

Algorithm:

- generate first times to negative and time to positive jump
- repeat 
    * if (time to last positive jump > time to last negative jump)
        # add last negative jump to path
        # if (stopping condition) exit loop
        # repeat
            % generate time to negative jump
            % if (time to last positive jump < time to last negative jump)
                @ exit loop
            % add negative jump to path
            % if (stopping condition) exit loop
        # if (stopping condition) exit loop
    * else 
        # add positive jump
        # if (stopping condition) exit loop
        # repeat 
            % generate time to positive jump
            % if (time to last positive jump > time to last negative jump)
                @ exit loop
            % add positive jump to path
            % if (stopping condition) exit loop
        # if (stopping condition) exit loop

Let me illustrate a couple of iterations to give a feeling of the algorithm. We simulate positive and negative jumps’ arrival:

Positive jump occurs later, therefore, we need to catch up with negative jumps. We simulate next negative arrival…

And next negative arrival…

One more…

And finally the negative jumps over positive and now we need to catch up with positve ones…

And so forth and so on. In the algorithm above stopping condition could be anythin, for instance, a maximum number of jumps is attained, the maximum number of iterations is attained, the maximum time span is attained, the path is ruined, etc. Below, I propose an implementation with stopping time that uses a maximum time span.

<span class="n">library</span><span class="p">(</span><span class="n">magrittr</span><span class="p">)</span><span class="w">
</span><span class="n">library</span><span class="p">(</span><span class="n">ggplot2</span><span class="p">)</span><span class="w">
</span><span class="n">library</span><span class="p">(</span><span class="n">skellam</span><span class="p">)</span><span class="w">

</span><span class="n">sim_p1</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">lambda_p</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">t</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">100</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
    
    </span><span class="c1"># utility function: get last element of a vector</span><span class="w">
    </span><span class="n">last</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="n">ifelse</span><span class="p">(</span><span class="nf">length</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="o">></span><span class="w"> </span><span class="m">0</span><span class="p">,</span><span class="w"> </span><span class="n">yes</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">x</span><span class="p">[</span><span class="nf">length</span><span class="p">(</span><span class="n">x</span><span class="p">)],</span><span class="w"> </span><span class="n">no</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w">
    
    </span><span class="c1"># initialize process</span><span class="w">
    </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">matrix</span><span class="p">(</span><span class="kc">NA</span><span class="p">,</span><span class="w"> </span><span class="n">nrow</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">ncol</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">2</span><span class="p">)</span><span class="w">
    </span><span class="n">colnames</span><span class="p">(</span><span class="n">path</span><span class="p">)</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="s2">"time"</span><span class="p">,</span><span class="w"> </span><span class="s2">"X"</span><span class="p">)</span><span class="w">
    </span><span class="n">path</span><span class="p">[</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="p">]</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="m">0</span><span class="p">,</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w">
    
    </span><span class="c1"># function for adding negative jump to a path</span><span class="w">
    </span><span class="n">add_jump_n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">()</span><span class="w"> </span><span class="p">{</span><span class="w">
        
        </span><span class="c1"># add a new time arrival to arrival times vector</span><span class="w">
        </span><span class="n">time_n</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="n">time_n</span><span class="p">,</span><span class="w"> </span><span class="n">current_time_n</span><span class="p">)</span><span class="w">
        
        </span><span class="c1"># add a negative jump to the path</span><span class="w">
        </span><span class="n">path</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
            </span><span class="n">path</span><span class="p">,</span><span class="w">
            </span><span class="nf">c</span><span class="p">(</span><span class="n">current_time_n</span><span class="p">,</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
        </span><span class="p">)</span><span class="w">
        
        </span><span class="n">path</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
            </span><span class="n">path</span><span class="p">,</span><span class="w">
            </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">-</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
        </span><span class="p">)</span><span class="w">
    </span><span class="p">}</span><span class="w">
    
    </span><span class="c1"># function for adding positive jump to a path</span><span class="w">
    </span><span class="n">add_jump_p</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">()</span><span class="w"> </span><span class="p">{</span><span class="w">
        
        </span><span class="c1"># add a new time arrival to arrival times vector</span><span class="w">
        </span><span class="n">time_p</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="n">time_p</span><span class="p">,</span><span class="w"> </span><span class="n">current_time_p</span><span class="p">)</span><span class="w">
        
        </span><span class="c1"># add a positive jump to the path</span><span class="w">
        </span><span class="n">path</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
            </span><span class="n">path</span><span class="p">,</span><span class="w">
            </span><span class="nf">c</span><span class="p">(</span><span class="n">current_time_p</span><span class="p">,</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
        </span><span class="p">)</span><span class="w">
        
        </span><span class="n">path</span><span class="w"> </span><span class="o"><<-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
            </span><span class="n">path</span><span class="p">,</span><span class="w">
            </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">+</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
        </span><span class="p">)</span><span class="w">
    </span><span class="p">}</span><span class="w">
    
    </span><span class="c1"># check whether the path is reached maximum time span</span><span class="w">
    </span><span class="n">is_max_time_span_attained</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">()</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">]</span><span class="w"> </span><span class="o">>=</span><span class="w"> </span><span class="n">t</span><span class="w">
    
    </span><span class="n">time_n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">numeric</span><span class="p">()</span><span class="w"> </span><span class="c1"># time of negative jumps arrivals</span><span class="w">
    </span><span class="n">time_p</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">numeric</span><span class="p">()</span><span class="w"> </span><span class="c1"># time of positive jumps arrivals</span><span class="w">
    
    </span><span class="n">current_time_n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rexp</span><span class="p">(</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_n</span><span class="p">)</span><span class="w"> </span><span class="c1"># current time arrival of the negative</span><span class="w">
    </span><span class="n">current_time_p</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rexp</span><span class="p">(</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_p</span><span class="p">)</span><span class="w"> </span><span class="c1"># current time arrival of the positive</span><span class="w">
    
    </span><span class="k">repeat</span><span class="p">{</span><span class="w">
        
        </span><span class="k">if</span><span class="p">(</span><span class="n">current_time_p</span><span class="w"> </span><span class="o">></span><span class="w"> </span><span class="n">current_time_n</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
            
            </span><span class="n">add_jump_n</span><span class="p">()</span><span class="w">
            
            </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            
            </span><span class="k">repeat</span><span class="w"> </span><span class="p">{</span><span class="w">
                
                </span><span class="n">current_time_n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">last</span><span class="p">(</span><span class="n">time_n</span><span class="p">)</span><span class="w"> </span><span class="o">+</span><span class="w"> </span><span class="n">rexp</span><span class="p">(</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_n</span><span class="p">)</span><span class="w">
                </span><span class="k">if</span><span class="p">(</span><span class="n">current_time_p</span><span class="w"> </span><span class="o"><</span><span class="w"> </span><span class="n">current_time_n</span><span class="p">)</span><span class="w"> </span><span class="k">break</span><span class="w">
                
                </span><span class="n">add_jump_n</span><span class="p">()</span><span class="w">
                
                </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            </span><span class="p">}</span><span class="w">
            
            </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            
            
        </span><span class="p">}</span><span class="w"> </span><span class="k">else</span><span class="w"> </span><span class="p">{</span><span class="w">
            
            </span><span class="n">add_jump_p</span><span class="p">()</span><span class="w">
            
            </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            
            </span><span class="k">repeat</span><span class="w"> </span><span class="p">{</span><span class="w">
                </span><span class="n">current_time_p</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">last</span><span class="p">(</span><span class="n">time_p</span><span class="p">)</span><span class="w"> </span><span class="o">+</span><span class="w"> </span><span class="n">rexp</span><span class="p">(</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_p</span><span class="p">)</span><span class="w">
                </span><span class="k">if</span><span class="p">(</span><span class="n">current_time_p</span><span class="w"> </span><span class="o">></span><span class="w"> </span><span class="n">current_time_n</span><span class="p">)</span><span class="w"> </span><span class="k">break</span><span class="w">
                
                </span><span class="n">add_jump_p</span><span class="p">()</span><span class="w">
                
                </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            </span><span class="p">}</span><span class="w">
            
            </span><span class="k">if</span><span class="p">(</span><span class="n">is_max_time_span_attained</span><span class="p">())</span><span class="w"> </span><span class="k">break</span><span class="w">
            
        </span><span class="p">}</span><span class="w">
    </span><span class="p">}</span><span class="w">
    
    </span><span class="c1"># dropping last step to be before t </span><span class="w">
    </span><span class="n">indices</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">path</span><span class="p">[,</span><span class="w"> </span><span class="m">1</span><span class="p">]</span><span class="w"> </span><span class="o"><=</span><span class="w"> </span><span class="n">t</span><span class="w">
    </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">indices</span><span class="p">,</span><span class="w"> </span><span class="p">,</span><span class="w"> </span><span class="n">drop</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="kc">FALSE</span><span class="p">]</span><span class="w">
    </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="n">path</span><span class="p">,</span><span class="w">
                  </span><span class="nf">c</span><span class="p">(</span><span class="n">t</span><span class="p">,</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]))</span><span class="w">
    
    </span><span class="n">rval</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">list</span><span class="p">(</span><span class="w">
        </span><span class="n">path</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">path</span><span class="p">,</span><span class="w">
        </span><span class="n">time_p</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">time_p</span><span class="p">,</span><span class="w">
        </span><span class="n">time_n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">time_n</span><span class="w">
    </span><span class="p">)</span><span class="w">
    
    </span><span class="nf">return</span><span class="p">(</span><span class="n">rval</span><span class="p">)</span><span class="w">
</span><span class="p">}</span><span class="w">
</span>

Let’s compare estimated expected value and variance of interarrival jumps with theoretical ones, for both positive and negative jumps. With default parameters $\lambda_1 = 1$ and $\lambda_2 = 1$, they all should be around one. This is confirmed in the code snipped below:

<span class="n">set.seed</span><span class="p">(</span><span class="m">2018</span><span class="p">)</span><span class="w">

</span><span class="n">p</span><span class="m">1</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sim_p1</span><span class="p">(</span><span class="n">t</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1000</span><span class="p">)</span><span class="w">
</span><span class="n">mean</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">1</span><span class="o">$</span><span class="n">time_p</span><span class="p">));</span><span class="w"> </span><span class="n">var</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">1</span><span class="o">$</span><span class="n">time_p</span><span class="p">))</span><span class="w">
</span><span class="c1"># [1] 0.9713893</span><span class="w">
</span><span class="c1"># [1] 0.9399382</span><span class="w">
</span><span class="n">mean</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">1</span><span class="o">$</span><span class="n">time_n</span><span class="p">));</span><span class="w"> </span><span class="n">var</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">1</span><span class="o">$</span><span class="n">time_n</span><span class="p">))</span><span class="w">
</span><span class="c1"># [1] 1.067944</span><span class="w">
</span><span class="c1"># [1] 1.121641</span><span class="w">
</span>

Method 2

The second method is a bit less cumbersome. We simulate separately the number of positive and negative jumps in the interval $(0, t)$ by Poisson distribution with respective rates. Then, we generate arrival times of jumps by the uniform distribution, which is then sorted. Finally, the full path should be built, that is adding positive and negative jumps in a loop. The idea is as before: we compare what kind of jump (negative vs positive) occure earlier, and add the earliest one. Note, that it is possible that several negative jumps occure before a positive one, and vice versa. Also it is possible that only positive or only negative jumps are left, therefore, we need to incorporate this in the loop.

The implementation is presetned below:

<span class="n">sim_p2</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">lambda_p</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda_n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">t</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">100</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
    
    </span><span class="c1"># simulate numbers of positive and negative jumps</span><span class="w">
    </span><span class="n">number_p_jumps</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rpois</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">lambda_p</span><span class="w"> </span><span class="o">*</span><span class="w"> </span><span class="n">t</span><span class="p">)</span><span class="w">
    </span><span class="n">number_n_jumps</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rpois</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">lambda</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">lambda_n</span><span class="w"> </span><span class="o">*</span><span class="w"> </span><span class="n">t</span><span class="p">)</span><span class="w">
    
    </span><span class="c1"># simulate the time of jumps' arrivals</span><span class="w">
    </span><span class="n">p_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">runif</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">number_p_jumps</span><span class="p">)</span><span class="w"> </span><span class="o">%>%</span><span class="w"> </span><span class="n">sort</span><span class="p">()</span><span class="w"> </span><span class="o">%>%</span><span class="w"> </span><span class="n">multiply_by</span><span class="p">(</span><span class="n">t</span><span class="p">)</span><span class="w">
    </span><span class="n">n_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">runif</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">number_n_jumps</span><span class="p">)</span><span class="w"> </span><span class="o">%>%</span><span class="w"> </span><span class="n">sort</span><span class="p">()</span><span class="w"> </span><span class="o">%>%</span><span class="w"> </span><span class="n">multiply_by</span><span class="p">(</span><span class="n">t</span><span class="p">)</span><span class="w">
    
    </span><span class="c1"># keep the time of jumps' arrivals in separate variables</span><span class="w">
    </span><span class="n">time_p</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">p_jumps_arrival</span><span class="w">
    </span><span class="n">time_n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">n_jumps_arrival</span><span class="w">
    
    </span><span class="c1"># initialize process</span><span class="w">
    </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">matrix</span><span class="p">(</span><span class="kc">NA</span><span class="p">,</span><span class="w"> </span><span class="n">nrow</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="n">ncol</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">2</span><span class="p">)</span><span class="w">
    </span><span class="n">colnames</span><span class="p">(</span><span class="n">path</span><span class="p">)</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="s2">"time"</span><span class="p">,</span><span class="w"> </span><span class="s2">"X"</span><span class="p">)</span><span class="w">
    </span><span class="n">path</span><span class="p">[</span><span class="m">1</span><span class="p">,</span><span class="w"> </span><span class="p">]</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">c</span><span class="p">(</span><span class="m">0</span><span class="p">,</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w">
    
    
    </span><span class="k">while</span><span class="p">(</span><span class="nf">length</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="w"> </span><span class="o">|</span><span class="w"> </span><span class="nf">length</span><span class="p">(</span><span class="n">n_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
        
        </span><span class="k">if</span><span class="p">(</span><span class="nf">length</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="w"> </span><span class="o">&</span><span class="w"> </span><span class="nf">length</span><span class="p">(</span><span class="n">n_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
            </span><span class="k">if</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">]</span><span class="w"> </span><span class="o"><</span><span class="w"> </span><span class="n">n_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">])</span><span class="w"> </span><span class="p">{</span><span class="w">
                
                </span><span class="c1"># add positive jump</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">+</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">p_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">p_jumps_arrival</span><span class="p">[</span><span class="m">-1</span><span class="p">]</span><span class="w">
                
            </span><span class="p">}</span><span class="w"> </span><span class="k">else</span><span class="w"> </span><span class="p">{</span><span class="w">
                
                </span><span class="c1"># add negative jump</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">n_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">-</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">n_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">n_jumps_arrival</span><span class="p">[</span><span class="m">-1</span><span class="p">]</span><span class="w">
                
            </span><span class="p">}</span><span class="w">
        </span><span class="p">}</span><span class="w"> </span><span class="k">else</span><span class="w"> </span><span class="p">{</span><span class="w">
            </span><span class="k">if</span><span class="p">(</span><span class="nf">length</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
                
                </span><span class="c1"># add positive jump</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">p_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">+</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">p_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">p_jumps_arrival</span><span class="p">[</span><span class="m">-1</span><span class="p">]</span><span class="w">
                
            </span><span class="p">}</span><span class="w">
            </span><span class="k">if</span><span class="p">(</span><span class="nf">length</span><span class="p">(</span><span class="n">n_jumps_arrival</span><span class="p">)</span><span class="w"> </span><span class="o">!=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
                
                </span><span class="c1"># add negative jump</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">n_jumps_arrival</span><span class="p">[</span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">])</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
                    </span><span class="n">path</span><span class="p">,</span><span class="w">
                    </span><span class="nf">c</span><span class="p">(</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">1</span><span class="p">],</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o">-</span><span class="w"> </span><span class="m">1</span><span class="p">)</span><span class="w">
                </span><span class="p">)</span><span class="w">
                
                </span><span class="n">n_jumps_arrival</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">n_jumps_arrival</span><span class="p">[</span><span class="m">-1</span><span class="p">]</span><span class="w">
            </span><span class="p">}</span><span class="w">
        </span><span class="p">}</span><span class="w">
    </span><span class="p">}</span><span class="w">
    
    </span><span class="c1"># add last step</span><span class="w">
    </span><span class="n">path</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="n">path</span><span class="p">,</span><span class="w">
                  </span><span class="nf">c</span><span class="p">(</span><span class="n">t</span><span class="p">,</span><span class="w"> </span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]))</span><span class="w">
    
    </span><span class="n">rval</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="nf">list</span><span class="p">(</span><span class="w">
        </span><span class="n">path</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">path</span><span class="p">,</span><span class="w">
        </span><span class="n">time_p</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">time_p</span><span class="p">,</span><span class="w">
        </span><span class="n">time_n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">time_n</span><span class="w">
    </span><span class="p">)</span><span class="w">
    
    </span><span class="nf">return</span><span class="p">(</span><span class="n">rval</span><span class="p">)</span><span class="w">
</span><span class="p">}</span><span class="w">
</span>

Again, we need to be assured that estimated expected values and variances for positive and negative jumps are close to one:

<span class="n">p</span><span class="m">2</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sim_p2</span><span class="p">(</span><span class="n">t</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1000</span><span class="p">)</span><span class="w">
</span><span class="n">mean</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">2</span><span class="o">$</span><span class="n">time_p</span><span class="p">));</span><span class="w"> </span><span class="n">var</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">2</span><span class="o">$</span><span class="n">time_p</span><span class="p">))</span><span class="w">
</span><span class="c1"># [1] 1.006591</span><span class="w">
</span><span class="c1"># [1] 1.054261</span><span class="w">
</span><span class="n">mean</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">2</span><span class="o">$</span><span class="n">time_n</span><span class="p">));</span><span class="w"> </span><span class="n">var</span><span class="p">(</span><span class="n">diff</span><span class="p">(</span><span class="n">p</span><span class="m">2</span><span class="o">$</span><span class="n">time_n</span><span class="p">))</span><span class="w">
</span><span class="c1"># [1] 0.9923275</span><span class="w">
</span><span class="c1"># [1] 0.8606385</span><span class="w">
</span>

Convergence

Finally, we want visually check if the estimatros are non-biased, and how fast they converge (i.e. which method has a smaller variance). We focus on the expected value of a path, as well as on the probability of the path to be below a certain value. For this we simulate a vast number (1000) of paths using each of methods. Then, we estimae the expected value as a function of the number of simulations by aggregating first x paths.

<span class="n">n</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="m">1000</span><span class="w">

</span><span class="n">paths1</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">replicate</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">expr</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">sim_p1</span><span class="p">(),</span><span class="w"> </span><span class="n">simplify</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="kc">FALSE</span><span class="p">)</span><span class="w">
</span><span class="n">paths2</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">replicate</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">expr</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">sim_p2</span><span class="p">(),</span><span class="w"> </span><span class="n">simplify</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="kc">FALSE</span><span class="p">)</span><span class="w">
</span>

First, we look at the expected value, which sould be zero given both lambdas equal one:

<span class="n">means1</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sapply</span><span class="p">(</span><span class="w">
    </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w">
    </span><span class="k">function</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
        </span><span class="n">mean</span><span class="p">(</span><span class="n">sapply</span><span class="p">(</span><span class="n">paths1</span><span class="p">[</span><span class="m">1</span><span class="o">:</span><span class="n">x</span><span class="p">],</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">y</span><span class="p">)</span><span class="w"> </span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]))</span><span class="w">
    </span><span class="p">}</span><span class="w">
</span><span class="p">)</span><span class="w">

</span><span class="n">means2</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sapply</span><span class="p">(</span><span class="w">
    </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w">
    </span><span class="k">function</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
        </span><span class="n">mean</span><span class="p">(</span><span class="n">sapply</span><span class="p">(</span><span class="n">paths2</span><span class="p">[</span><span class="m">1</span><span class="o">:</span><span class="n">x</span><span class="p">],</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">y</span><span class="p">)</span><span class="w"> </span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]))</span><span class="w">
    </span><span class="p">}</span><span class="w">
</span><span class="p">)</span><span class="w">

</span><span class="n">means</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
    </span><span class="n">data.frame</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">mean</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">means1</span><span class="p">,</span><span class="w"> </span><span class="n">method</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="s2">"1"</span><span class="p">),</span><span class="w">
    </span><span class="n">data.frame</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">mean</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">means2</span><span class="p">,</span><span class="w"> </span><span class="n">method</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="s2">"2"</span><span class="p">)</span><span class="w">
</span><span class="p">)</span><span class="w">

</span><span class="n">ggplot</span><span class="p">(</span><span class="n">means</span><span class="p">)</span><span class="w"> </span><span class="o">+</span><span class="w">
    </span><span class="n">geom_line</span><span class="p">(</span><span class="n">aes</span><span class="p">(</span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">mean</span><span class="p">,</span><span class="w"> </span><span class="n">color</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">method</span><span class="p">))</span><span class="w"> </span><span class="o">+</span><span class="w">
    </span><span class="n">geom_hline</span><span class="p">(</span><span class="n">yintercept</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">0</span><span class="p">)</span><span class="w"> </span><span class="o">+</span><span class="w"> 
    </span><span class="n">theme_bw</span><span class="p">()</span><span class="w"> </span><span class="o">+</span><span class="w"> 
    </span><span class="n">theme</span><span class="p">(</span><span class="n">text</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">element_text</span><span class="p">(</span><span class="n">size</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">24</span><span class="p">))</span><span class="w">
</span>

From the plot it seems that both methods converge with the same speed to the correct value. Huray! How about probabilities?

For probabilites we used a threashold of 10 (arbitrary choosen), and a true value calculated by pskellam package. For default argument t = 10, the distribution of $X(10) \sim Skellam(\lambda_1 \cdot t, \lambda_2 \cdot t)$.

<span class="n">probs1</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sapply</span><span class="p">(</span><span class="w">
    </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w">
    </span><span class="k">function</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
        </span><span class="n">mean</span><span class="p">(</span><span class="n">sapply</span><span class="p">(</span><span class="n">paths1</span><span class="p">[</span><span class="m">1</span><span class="o">:</span><span class="n">x</span><span class="p">],</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">y</span><span class="p">)</span><span class="w"> </span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o"><=</span><span class="w"> </span><span class="m">10</span><span class="p">))</span><span class="w">
    </span><span class="p">}</span><span class="w">
</span><span class="p">)</span><span class="w">

</span><span class="n">probs2</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">sapply</span><span class="p">(</span><span class="w">
    </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w">
    </span><span class="k">function</span><span class="p">(</span><span class="n">x</span><span class="p">)</span><span class="w"> </span><span class="p">{</span><span class="w">
        </span><span class="n">mean</span><span class="p">(</span><span class="n">sapply</span><span class="p">(</span><span class="n">paths2</span><span class="p">[</span><span class="m">1</span><span class="o">:</span><span class="n">x</span><span class="p">],</span><span class="w"> </span><span class="k">function</span><span class="p">(</span><span class="n">y</span><span class="p">)</span><span class="w"> </span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">[</span><span class="n">nrow</span><span class="p">(</span><span class="n">y</span><span class="o">$</span><span class="n">path</span><span class="p">),</span><span class="w"> </span><span class="m">2</span><span class="p">]</span><span class="w"> </span><span class="o"><=</span><span class="w"> </span><span class="m">10</span><span class="p">))</span><span class="w">
    </span><span class="p">}</span><span class="w">
</span><span class="p">)</span><span class="w">

</span><span class="n">probs</span><span class="w"> </span><span class="o"><-</span><span class="w"> </span><span class="n">rbind</span><span class="p">(</span><span class="w">
    </span><span class="n">data.frame</span><span class="p">(</span><span class="n">n</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="m">1</span><span class="o">:</span><span class="n">n</span><span class="p">,</span><span class="w"> </span><span class="n">prob</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="n">probs1</span><span class="p">,</span><span class="w"> </span><span class="n">method</span><span class="w"> </span><span class="o">=</span><span class="w"> </span><span class="s2">"1"</span><span class="p">),</span><span class="w">
    </span><span class="n">data.frame</span><span class="...

To leave a comment for the author, please follow the link and comment on their blog: Iegor Rudnytskyi.

R-bloggers.com offers daily e-mail updates about R news and tutorials about learning R and many other topics. Click here if you're looking to post or find an R/data-science job.
Want to share your content on R-bloggers? click here if you have a blog, or here if you don't.

Never miss an update!
Subscribe to R-bloggers to receive
e-mails with the latest R posts.
(You will not see this message again.)

Click here to close (This popup will not appear again)