PyHawkes/GenuineMultivariateHawkesProcessFitter.py at master · thomasj02/PyHawkes · GitHub

1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
137
138
139
140
__author__ = 'tjohnson'
import ImmigrationDescendantParameters
import random
import scipy.optimize

class GenuineMultivariateHawkesProcessFitter:
    def __init__(self,hawkesProcess):
        self.hawkesProcess=hawkesProcess

        totalNumParams=0
        parameterBounds=[]

        numComponents=hawkesProcess.numComponents
        self.immigrationDescendantNumParams=ImmigrationDescendantParameters.ImmigrationDescendantParameters.getNumParameters(numComponents)
        totalNumParams+=self.immigrationDescendantNumParams
        parameterBounds.extend(ImmigrationDescendantParameters.ImmigrationDescendantParameters.getParameterBounds(numComponents))

        self.decayFunctionsNumParams=[]
        for decayFunction in self.__getDecayFunctions():
            thisDecayFunctionNumParams=decayFunction.getNumParameters()
            self.decayFunctionsNumParams.append(thisDecayFunctionNumParams)
            totalNumParams+=thisDecayFunctionNumParams
            parameterBounds.extend(decayFunction.getParameterBounds())

        self.markDistributionsNumParams=[]
        for markDistribution in self.__getMarkDistributions():
            thisMarkDistributionNumParams=markDistribution.getNumParameters()
            self.markDistributionsNumParams.append(thisMarkDistributionNumParams)
            totalNumParams+=thisMarkDistributionNumParams
            parameterBounds.extend(markDistribution.getParameterBounds())

        self.numParams=totalNumParams
        self.parameterBounds=parameterBounds

    def __getMarkDistributions(self):
        retval=[]

        lastMarkDistribution=None
        for markDistribution in self.hawkesProcess.markDistributions:
            if markDistribution is lastMarkDistribution:
                continue
            retval.append(markDistribution)
            lastMarkDistribution=markDistribution

        return retval

    def __getDecayFunctions(self):
        retval=[]

        lastDecayFunction=None
        for decayFunction in self.hawkesProcess.decayFunctions:
            if decayFunction is lastDecayFunction:
                continue
            retval.append(decayFunction)
            lastDecayFunction=decayFunction

        return retval

    def getInitialRandomVector(self):
        vector=[]
        for lowerBound,upperBound in self.parameterBounds:
            if lowerBound is None:
                lowerBound=0.0
            randomValue=random.uniform(lowerBound,lowerBound+1.0)
            vector.append(randomValue)

        return vector

    def setParameterValues(self,params):
        parameterStartIndex=0
        self.hawkesProcess.immigrationDescendantParameters.setParameters(params[parameterStartIndex:self.immigrationDescendantNumParams])
        parameterStartIndex+=self.immigrationDescendantNumParams

        for decayFunction,decayFunctionNumParams in zip(self.__getDecayFunctions(),self.decayFunctionsNumParams):
            decayFunctionParams=params[parameterStartIndex:parameterStartIndex+decayFunctionNumParams]
            parameterStartIndex+=decayFunctionNumParams

            decayFunction.setParams(decayFunctionParams)


        for markDistribution,markDistributionNumParams in zip(self.__getMarkDistributions(),self.markDistributionsNumParams):
            markDistributionParams = params[parameterStartIndex:parameterStartIndex + markDistributionNumParams]
            parameterStartIndex+=markDistributionNumParams

            markDistribution.setParams(markDistributionParams)

    def __getNegativeLogLikelihoodWithParameters(self,parameters,timeComponentMarkTriples):
        print "Getting log likelihood with parameters",parameters
        self.setParameterValues(parameters)
        negativeLogLikelihood=-self.hawkesProcess.getLogLikelihood(timeComponentMarkTriples)
        print "Got",negativeLogLikelihood
        return negativeLogLikelihood

    def maximizeLikelihood(self,timeComponentMarkTriples,initialGuess):
        self.setParameterValues(initialGuess)
        negativeLogLikelihoodFunction=lambda x: self.__getNegativeLogLikelihoodWithParameters(x,timeComponentMarkTriples)
        print "Initial Likelihood: %s" % negativeLogLikelihoodFunction(initialGuess)

        bestParams=scipy.optimize.fmin_l_bfgs_b(negativeLogLikelihoodFunction,x0=initialGuess,approx_grad=True,bounds=self.parameterBounds,iprint=1)
        print "BestParams:",bestParams
        print "Best Log likelihood:",negativeLogLikelihoodFunction(bestParams)


if __name__=="__main__":
    import DecayFunctions
    import MarkDistributions
    import GenuineMultivariateHawkesProcess
    import ImmigrationDescendantParameters
    import numpy as np


    q=np.matrix("[0.61 0.16;0.60 0.06]")
    nu=np.array([0.021,0.029])
    alpha=0.015
    rho1=5.6
    rho2=7.2
    mu1=3.6
    mu2=4.2
    phi1=0.47
    phi2=1.1
    psi1=0.22
    psi2=0.0

    immigrationDescendantParameters=\
        ImmigrationDescendantParameters.ImmigrationDescendantParameters(2,[0.021,0.029,0.61,0.16,0.6,0.06])
    print immigrationDescendantParameters.nu
    print immigrationDescendantParameters.q
    decayFunction=DecayFunctions.ExponentialDecayFunction([alpha])
    markDistribution1=MarkDistributions.ParetoMarkDistribution([mu1,rho1,phi1,psi1,0.0])
    markDistribution2=MarkDistributions.ParetoMarkDistribution([mu2,rho2,phi2,psi2,0.0])

    hawkesProcess=GenuineMultivariateHawkesProcess.GenuineMultivariateHawkesProcess(
        immigrationDescendantParameters,
        [decayFunction,decayFunction],
        [markDistribution1,markDistribution2])

    fitter=GenuineMultivariateHawkesProcessFitter(hawkesProcess)
    print fitter.numParams
    print fitter.parameterBounds
    print fitter.getInitialRandomVector()