Aide pour finaliser ce code
Hello Ablaye, I didn't find any function xi ,yi. This is an english mailing list. You will maximize your chance of getting help by sticking to english. I would also recommend that you comment your code in English, it will be easier to share. But on this one it is your call. PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code. Best regards, Jeremie Juste
On Friday, 9 Oct 2020 at 13:31, Ablaye Ngalaba wrote:
Hello.
Here is my R code. I used the functional data . Now I need to use the
functional data by applying the kernels instead of the xi, yi functions.
Bonjour.
Voici mon code en R . J'ai utiliser les donn?es fonctionnelles . Maintenant
j'ai besoin d'utiliser les donn?es fonctionnelles en appliquant les noyaux
? la place des fontions xi, yi
library(MASS)
CentrageV<-function(X,Ms,n){
# cette fonction centre les donn?es de X
X1=X*0
for (i in 1:n){
X1[i,]=X[i,]-Ms
}
return(X1)
}
# Fonction N?2
SqrtMatrice0<-function(M) {
# Cette fonction nous permet de calculer la racine carr?e d'une matrice
# en utilisant la d?composition M=PDQ o? Q est l'inverse de P
# ici les valeurs propres n?gatives sont remplac?es par zero
a=eigen(M,TRUE)
b=a$values
b[b<0]=0
c=a$vectors
d=diag(sqrt(b))
b=solve(c)
a=c%*%d%*%b
return(a)
}
# d?claration des parametres
m1=0.01 # valeur de alpha (risque de 1%)
m2=0.05 # valeur de alpha (risque de 5%)
m3=0.1 # valeur de alpha (risque de 10%)
nbrefoissim=100 # nbrefois que le programme tourne
p=2 #dimension de la variable X
q=3 #dimension de la variable Y
R=c(2,3,2);# Nbre de partition de chaque composante de la variable Y
if(length(R) != q) stop("La taille de R doit ?tre ?gale ? q")
n=25 # Taille echantillon
N=c(25,50,100,200,300,400,500,1000) #differentes tailles de l'?chantillion
N=c(25,50) #differentes tailles de l'?chantillion
K=0
MV=matrix(0,nr=2,nc=4)
dimnames(MV)[[2]]=c("n","r1%","r5%","r10%")
#Debut du programme
for (n in N){
l1=0 # initialisation de la valeur permettant calculer le niveau de test ?
1%
l2=0 # initialisation de la valeur permettant calculer le niveau de test ?
5%
l3=0 # initialisation de la valeur permettant calculer le niveau de test ?
10%
# Cr?ation d'une liste n11 qui contient les tailles des differents groupes
n11=list()
for (i in 1:q){
n11[[i]]=rep(as.integer(n/R[i]),R[i])
n11[[i]][R[i]]=n-((R[i]-1)*n11[[i]][1])
}
# Cr?ation des listes P11 et P12 qui contient les probabilit?s et
# les inverses des probabilites empiriques des differents groupes
respectivement
P11=list()
P12=list()
for (i in 1:q){
P11[[i]]=n11[[i]]/n
P12[[i]]=n/n11[[i]]
}
# cr?ation d'une liste contenant les matrices W
W=list()
for (i in 1:q){
w=matrix(0,n,R[i])
w[1:n11[[i]][1],1]=1
for (j in 2:R[i]){
s1=sum(n11[[i]][1:(j-1)])
w[(1+s1):(s1+n11[[i]][j]),j]=1
}
W[[i]]=w
}
for (i1 in 1:nbrefoissim){
# g?neration des donn?es
VA1=mvrnorm(n,rep(0,(p+q)),diag((p+q)))
X=VA1[,1:p]
Y=VA1[,(p+1):(p+q)]
# Calcul de Xbar
Xbar=colMeans(X)
# Calcul des Xjh bar
Xjhbar=list()
for (i in 1:q){
w=matrix(0,R[i],p)
for (j in 1:R[i]){
w[j,]=colSums(W[[i]][,j]*X)/n11[[i]][j]
}
Xjhbar[[i]]=w
}
#calcul des TO jh
TO.jh=list()
for (i in 1:q){
w=Xjhbar[[i]]
to=w*0
for (j in 1:R[i]){
to[j,]=w[j,]-Xbar
}
TO.jh[[i]]=to
}
#calcul des Lamda J
Lamda=matrix(0,p,p)
for (i in 1:q){
to=TO.jh[[i]]
w=matrix(0,p,p)
for (j in 1:R[i]){
w=w+(P11[[i]][j]*(to[j,]%*%t(to[j,])))
}
Lamda=Lamda+w
}
tr1=n*sum(diag(Lamda))
# Calcul de Gamma
GGamma=matrix(0,p*sum(R),p*sum(R))
PGamma=kronecker(diag(P11[[1]]),diag(p))
Ifin=p*R[1]
GGamma[1:Ifin,1:Ifin]=PGamma
for (i in 2:q){
PGamma=kronecker(diag(P11[[i]]),diag(p))
Idebut=((p*sum(R[1:(i-1)]))+1)
Ifin=(p*sum(R[1:i]))
GGamma[Idebut:Ifin,Idebut:Ifin]=PGamma
}
#Calcul de Sigma
# Calcul de Vn
X1=CentrageV(X,Xbar,n)
Vn=t(X1)%*%X1/n
## Construction de Sigma
GSigma=matrix(0,p*sum(R),p*sum(R))
for (i in 1:q ){
for (j in 1:R[i] ){
for (k in 1:q ){
for (l in 1:R[k]){
Xij=CentrageV(X,Xjhbar[[i]][j,],n)
Xkl=CentrageV(X,Xjhbar[[k]][l,],n)
Vijkl=t(W[[i]][j]*Xij)%*%(W[[k]][l]*Xkl)/n
Vij=t(W[[i]][j]*Xij)%*%Xij/n11[[i]][j]
Vkl=t(W[[k]][l]*Xkl)%*%Xkl/n11[[k]][l]
if (i==1) Idebut=((j-1)*p)+1 else
Idebut=((sum(R[1:(i-1)])+j-1)*p)+1
if (i==1) Idebut1=(j*p) else Idebut1=((sum(R[1:(i-1)])+j)*p)
if (k==1) Ifin=((l-1)*p)+1 else Ifin=((sum(R[1:(k-1)])+l-1)*p)+1
if (k==1) Ifin1=(l*p) else Ifin1=((sum(R[1:(k-1)])+l)*p)
GSigma[Idebut:Idebut1,Ifin:Ifin1]=Vijkl+(P11[[i]][j]*P11[[k]][l]*(Vn-Vij-Vkl))
}
}
}
}
# D?terminations des valeurs propres de sigma epsilon
pa=SqrtMatrice0(GSigma)
mq= pa %*% GGamma %*% pa
u=Re(eigen(mq)$values)
# d?termination de d?gr? de libert? et valeur c not? va
dl=(sum(u)^2)/(sum(u^2))
va=(sum(u^2))/(sum(u))
pc=1-pchisq(tr1/va, df= dl)
# Test de la valeur obtenue
if (pc>m1) d1=0 else d1=1
if (pc>m2) d2=0 else d2=1
if (pc>m3) d3=0 else d3=1
l1=l1+d1
l2=l2+d2
l3=l3+d3
}
K=K+1
MV[K,1]=n
MV[K,2]=l1/nbrefoissim
MV[K,3]=l2/nbrefoissim
MV[K,4]=l3/nbrefoissim
}
[[alternative HTML version deleted]]
______________________________________________ R-help at r-project.org mailing list -- To UNSUBSCRIBE and more, see https://stat.ethz.ch/mailman/listinfo/r-help PLEASE do read the posting guide http://www.R-project.org/posting-guide.html and provide commented, minimal, self-contained, reproducible code.