From 2e89ba20b5c9e067ef8eaf80e6071c1f7906510d Mon Sep 17 00:00:00 2001 From: Paul-Corbalan <58653590+Paul-Corbalan@users.noreply.github.com> Date: Tue, 5 Apr 2022 10:50:23 +0200 Subject: [PATCH] Local score work --- Comparaison_of_methods.rmd | 37 +++++++++++++++++++++---------------- 1 file changed, 21 insertions(+), 16 deletions(-) diff --git a/Comparaison_of_methods.rmd b/Comparaison_of_methods.rmd index 86f8c92..bd428c7 100644 --- a/Comparaison_of_methods.rmd +++ b/Comparaison_of_methods.rmd @@ -123,12 +123,12 @@ PValue <- function(Emp,ppH1, T, tau){ ### 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. ```{r} -NbSeqH0=5000 +NbSeqH0=1000 NbSeqH1=NbSeqH0 DataH0=vector("list") DataH1=vector("list") lambda0=3 -lambda1=4 +lambda1=5 T=10 tau=1 @@ -164,7 +164,7 @@ We compute the p-value associated to all 5 sequences, and stock them in a vector ```{r} #We start by computing the empirical distribution for lambda0 Emp=EmpDistrib(lambda0,n_sample,T,tau) -scan=c() + pvalue=c() index_scan=c() @@ -220,15 +220,23 @@ ppH0 = PoissonProcess(lambda0,10^4) n1 = length(ppH0) tbe0 = ppH0[2:n1]-ppH0[1:n1-1] 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) #print(summary(x)) -P_X = table(factor(x, levels = min(x):max(x))) -P_X = P_X/sum(table(x)) +min_X = min(xp) +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(P_X) +#print(P_X) ``` ### Calcul du local score @@ -237,19 +245,16 @@ library("localScore") library(Rcpp) E = 10 pvalue=c() +X=c() 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) - min_X = min(X) - P_X = table(factor(X, levels = min_X:max_X))/length(X) # supprimer pour les séquences de MC + result = daudin(localScore = LS, score_probabilities = P_X, sequence_length = length(x), sequence_min = min_X, sequence_max = max_X) - LS=localScoreC(X)$localScore[1] - - result = daudin(localScore = LS, score_probabilities = P_X, sequence_length = length(X), sequence_min = min_X, sequence_max =max_X) - - pvalue=c(pvalue,result) + pvalue = c(pvalue,result) } LS_H0=data.frame(num=1:NbSeqH0, pvalue_scan=pvalue, class=(pvalue<0.05)) LS_H0