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
|
#syntaxe de la log vraisemblance tirée à partir de la loi de pareto
pareto.like<-function(C,alpha,y){
y=as.matrix(y);
n<-nrow(y);
logl<-log(C)-log(alpha)+((alpha+1)/C)*sum(log(y))
return(logl)
}
#max de vraisemblance
#optim(starting values, function, data, method)
#optimisation de la Loi de pareto
#Pour cela on va simuler les pareto par méthode de l'inversion (rpareto ne marche pas dans notre cas)
invfdrpareto=function(y,C,a)
{x=(C/(1-y))^(1/a)
x}
rpareto=function(n,C,a)
{U=runif(n)
Y=sapply(U,invfdrpareto,C,a)
Y}
y=rpareto(1000,1,1000)
#log vraisemblance d'une pareto
logy=pareto.like(1,1000,y)
a=optim(1,pareto.like,y=y,method="BFGS") |
Partager