paul-corbalan
/
Scan-Statistics-Project-4Y-INSA
1
0
Fork 0
Code Issues Pull Requests Packages Projects Releases Wiki Activity

Local score work

This commit is contained in:
Paul-Corbalan 2022-04-05 10:50:23 +02:00
parent 7e8eb49d95
commit 2e89ba20b5
1 changed files with 21 additions and 16 deletions
Show Stats Download Patch File Download Diff File Expand all files Collapse all files

37
Comparaison_of_methods.rmd
View File

@ -123,12 +123,12 @@ PValue <- function(Emp,ppH1, T, tau){
### 2.2. Simulation under $\mathcal{H}_0$ and computation of p-values ### 2.2. Simulation under $\mathcal{H}_0$ and computation of p-values
On simule des séquences sous $\mathcal{H}_0$, que l'on stocke. On calcule la valeur de la scan stat et de la p-value, que l'on stocke aussi. On a une séquence de p-valeur des scans et une séquence de score local. On simule des séquences sous $\mathcal{H}_0$, que l'on stocke. On calcule la valeur de la scan stat et de la p-value, que l'on stocke aussi. On a une séquence de p-valeur des scans et une séquence de score local.
```{r} ```{r}
NbSeqH0=5000 NbSeqH0=1000
NbSeqH1=NbSeqH0 NbSeqH1=NbSeqH0
DataH0=vector("list") DataH0=vector("list")
DataH1=vector("list") DataH1=vector("list")
lambda0=3 lambda0=3
lambda1=4 lambda1=5
T=10 T=10
tau=1 tau=1
@ -164,7 +164,7 @@ We compute the p-value associated to all 5 sequences, and stock them in a vector
```{r} ```{r}
#We start by computing the empirical distribution for lambda0 #We start by computing the empirical distribution for lambda0
Emp=EmpDistrib(lambda0,n_sample,T,tau) Emp=EmpDistrib(lambda0,n_sample,T,tau)
scan=c()
pvalue=c() pvalue=c()
index_scan=c() index_scan=c()
@ -220,15 +220,23 @@ ppH0 = PoissonProcess(lambda0,10^4)
n1 = length(ppH0) n1 = length(ppH0)
tbe0 = ppH0[2:n1]-ppH0[1:n1-1] tbe0 = ppH0[2:n1]-ppH0[1:n1-1]
print(ks.test(tbe0,'exp')) print(ks.test(tbe0,'exp'))
x = floor(E*(log(lambda1/lambda0)+(lambda0-lambda1)*tbe0)) # ne pas mettre le floor ni le E (certes égale à 1) xp = floor(E*(log(lambda1/lambda0)+(lambda0-lambda1)*tbe0)) # ne pas mettre le floor ni le E (certes égale à 1)
#hist(x) #hist(x)
#print(summary(x)) #print(summary(x))
P_X = table(factor(x, levels = min(x):max(x))) min_X = min(xp)
P_X = P_X/sum(table(x)) max_X = max(xp)
vect.score = min_X:max_X
P_X = table(factor(xp, levels = min(xp):max(xp)))
P_X = P_X/sum(table(xp))
Mean_xp = sum(vect.score*P_X)
print(Mean_xp)
#print(dist.theo.scores) # Mettre à jour avec Elisa #print(dist.theo.scores) # Mettre à jour avec Elisa
print(P_X) #print(P_X)
``` ```
### Calcul du local score ### Calcul du local score
@ -237,19 +245,16 @@ library("localScore")
library(Rcpp) library(Rcpp)
E = 10 E = 10
pvalue=c() pvalue=c()
X=c()
for (i in 1:NbSeqH0){ for (i in 1:NbSeqH0){
X = floor(E*log(dexp(tbe0[[i]], rate = lambda1)/dexp(tbe0[[i]], rate = lambda0))) x = floor(E*log(dexp(tbe0[[i]], rate = lambda1)/dexp(tbe0[[i]], rate = lambda0)))
X=c(X,x)
LS=localScoreC(x)$localScore[1]
max_X = max(X) result = daudin(localScore = LS, score_probabilities = P_X, sequence_length = length(x), sequence_min = min_X, sequence_max = max_X)
min_X = min(X)
P_X = table(factor(X, levels = min_X:max_X))/length(X) # supprimer pour les séquences de MC
LS=localScoreC(X)$localScore[1] pvalue = c(pvalue,result)
result = daudin(localScore = LS, score_probabilities = P_X, sequence_length = length(X), sequence_min = min_X, sequence_max =max_X)
pvalue=c(pvalue,result)
} }
LS_H0=data.frame(num=1:NbSeqH0, pvalue_scan=pvalue, class=(pvalue<0.05)) LS_H0=data.frame(num=1:NbSeqH0, pvalue_scan=pvalue, class=(pvalue<0.05))
LS_H0 LS_H0