Hallo everyone!

I'm a beginner in PERL and have been trying to implement in PERL a code that would simulate a uni-dimensional Hawkess process. I have found several algorithms to do so and decided to go with the one suggested by Ogata (1981) which appears on slides 21/22 of the following presentation:

http://fiquant.mas.ecp.fr/ioane_files/HawkesCourseSlides.pdf

the code seems to work fine but when I try to use MLE estimation to get the values of the parameters used to generate the sequence, the base intensity is always way off while the jump and decay parameters are estimated accurately... why could this be happening? is there a mistake in my code? thanks in advance for how ever will take their time to help me out! regards

kris

I use R to perform the MLE and the following is my code in PERL: 

    use strict;
    use warnings;
    use diagnostics;
    
    use Data::Dumper;
    use List::Util qw(sum);
    use Math::Random::MT qw(rand);
    
    
    open (F1,'>>event_times_new2.csv');
    open (F2,'>>lambda_star_new2.csv');
    
    my @event_times;
    my @net_event_times;
    my @immigrants;
    my @new_immigrants;
    my @s;
    my @d;
    my @u;
    my @v;
    my $lambda0;
    my $alpha;
    my $beta;
    my $result;
    my $sum;
    my $lambda_star;
    my $time;
    
    
    my $counterfails=0;
    my $num=1;
    
    
    $time=8000000;
    
    $lambda0=1.2;
    $alpha=0.8;
    $beta=0.9;
    
    #subroutine that computes intensity function
    
    sub lambda
    {
    my $t = $_[0];
    $sum=0;
    for (my $i=0; $i<@event_times; $i++)
    {
    
    $net_event_times[$i]=($t-$event_times[$i]);
    
    $immigrants[$i]=$alpha*exp(((-1)*$beta)*$net_event_times[$i]);
    
    }
    
    $sum= sum @immigrants;
    
    if ($immigrants[-1]==$alpha){
    $sum+=(-$alpha);
    }
    $result=$lambda0+$sum;
    
    return $result;
    }
    
    
    #ALGO it self
    
    
    for (my $i=0; $i<2; $i++)
    {
    if ($i==0){
    $lambda_star=$lambda0;
    $s[$i]=0;
    $d[$i]=0;
    
    $u[$i]=rand($num);
    $s[$i]=(-1)*(1/$lambda_star)*log($u[$i]);
    if ($s[$i]<$time){
    $event_times[$i]=$s[$i];
    
    print F1 "$event_times[$i]\n";
    print F2 "$lambda_star\n";
    }else{
    die
    }
    }else
    {
    my $valuelambda=lambda($event_times[(($i)-1)]);
    $lambda_star=$valuelambda+$alpha;
    
    START:
    $u[$i]=rand($num);
    $s[$i]=$s[(($i)-1)]+(-1)*(1/$lambda_star)*log($u[$i]);
    if ($s[$i]>$time){
    die
    }
    
    $d[$i]=rand($num) ;
    my $rejectionlambda=lambda($s[$i]);
    if ($d[$i]<($rejectionlambda/$lambda_star)){
    $event_times[$i]=$s[$i];
    
    print F1 "$event_times[$i]\n";
    print F2 "$lambda_star\n";
    }else{
    
    $counterfails+=1;
    $lambda_star=$rejectionlambda;
    goto START;
    }
    }
    }
    print F1 "The number of event-times that failed the A/R procedure 
+is $counterfails\n";
    
    close F1;
    close F2;

[download]

In reply to Simulating uni-dimensional Hawkes processes by glrm_master

Title:
Use:  <p> text here (a paragraph) </p>
and:  <code> code here </code>
to format your post, it's "PerlMonks-approved HTML":


  • Posts are HTML formatted. Put <p> </p> tags around your paragraphs. Put <code> </code> tags around your code and data!
  • Titles consisting of a single word are discouraged, and in most cases are disallowed outright.
  • Read Where should I post X? if you're not absolutely sure you're posting in the right place.
  • Please read these before you post! —
  • Posts may use any of the Perl Monks Approved HTML tags:
    a, abbr, b, big, blockquote, br, caption, center, col, colgroup, dd, del, details, div, dl, dt, em, font, h1, h2, h3, h4, h5, h6, hr, i, ins, li, ol, p, pre, readmore, small, span, spoiler, strike, strong, sub, summary, sup, table, tbody, td, tfoot, th, thead, tr, tt, u, ul, wbr
  • You may need to use entities for some characters, as follows. (Exception: Within code tags, you can put the characters literally.)
            For:     Use:
    & &amp;
    < &lt;
    > &gt;
    [ &#91;
    ] &#93;
  • Link using PerlMonks shortcuts! What shortcuts can I use for linking?
  • See Writeup Formatting Tips and other pages linked from there for more info.