Update of distribtion via MC

This commit is contained in:
Paul-Corbalan 2022-04-03 18:31:36 +02:00
parent 425077cf27
commit fb66d1307a
1 changed files with 22 additions and 24 deletions

View File

@ -189,42 +189,40 @@ LS_H0
``` ```
## A reformater ## A reformater
```{r}
# distribtion des scores via MC # distribtion des scores via MC
# Nb seq. pp -> Nb seq. tbe -> dist. tbe (vérif) + Nb seq. Scores -> distr scores # Nb seq. pp -> Nb seq. tbe -> dist. tbe (vérif) + Nb seq. Scores -> distr scores
PoissonProcess <- function(lambda,T) { A = 1/(lambda0-lambda1)
return(sort(runif(rpois(1,lambda*T),0,T))) B = A*log(lambda1/lambda0)
}
lambda0=2 ppH1 = PoissonProcess(lambda1,T)
lambda1=3 n1 = length(ppH1)
Ti=100000 tbe1 = ppH1[2:n1]-ppH1[1:n1-1]
pp1=PoissonProcess(lambda1,Ti) print(tbe1)
print(pp1) print(ks.test(tbe1,'exp'))
n1=length(pp1) x = log(lambda1/lambda0)+(lambda0-lambda1)*tbe1 # ne pas mettre le floor ni le E (certes égale à 1)
tbe1=pp1[2:n1]-pp1[1:n1-1]
tbe1
ks.test(tbe1,'exp',lambda)
x=log(lambda1/lambda0)+(lambda0-lambda1)*tbe1 # ne pas mettre le floor ni le E (certes égale à 1)
hist(x) hist(x)
summary(x) print(summary(x))
# Calcul du maximum des scores # Calcul du maximum des scores
E=1 E = 1
# THEO à faire !!! max.s=log(lambda1/lambda0) # THEO à faire !!! max.s = log(lambda1/lambda0)
maxXk = floor(E*(log(lambda1/lambda0))) maxXk = floor(E*(log(lambda1/lambda0)))
maxXk maxXk
while (maxXk < 3) { while (maxXk < 3) {
E = E+1 E = E+1
maxXk = floor(E*(log(lambda1/lambda0))) maxXk = floor(E*(log(lambda1/lambda0)))
} }
E print(E)
x=floor(E*(log(lambda1/lambda0)+(lambda0-lambda1)*tbe1)) x = floor(E*(log(lambda1/lambda0)+(lambda0-lambda1)*tbe1))
dist.emp.scores=table(x)/sum(table(x)) dist.emp.scores = table(x)/sum(table(x))
dist.emp.scores dist.emp.scores
hist(x) hist(x)
range(x) print(range(x))
x.verif=seq(range(x)[1],range(x)[2],1) x.verif = seq(range(x)[1],range(x)[2],1)
dist.theo.scores=lambda0*exp(-lambda0*(A*x.verif-B)) dist.theo.scores = lambda0*exp(-lambda0*(A*x.verif-B))
dist.theo.scores print(dist.theo.scores)
dist.emp.scores print(dist.emp.scores)
```