R/SEML.R

In musto101/wilcox_R: What the package does (short line)

SEML <-
function(Z,thetaMLE){
	n=nrow(Z)
	p=ncol(Z)
	ps=p*(p+1)/2
   	q=ps
	zbar=MeanCov(Z)$zbar
	S=MeanCov(Z)$S
	Dup=Dp(p)
	InvS=solve(S)
	W=0.5*t(Dup)%*%(InvS%x%InvS)%*%Dup
	OmegaInf=solve(W)  # only about sigma, not mu
	
	
	# Sandwich-type Omega
	S12=matrix(0,p,ps)
	S22=matrix(0,ps,ps)
	
	for (i in 1:n){
		zi0=Z[i,]-zbar
		difi=zi0%*%t(zi0)-S
   	   	vdifi=vech(difi)   		
		S12=S12+zi0%*%t(vdifi)	
		S22=S22+vdifi%*%t(vdifi)		
	}
	
	OmegaSW=S22/n # only about sigma, not mu
	
   	SEinf=getSE(thetaMLE,OmegaInf,n)
   	SEsw=getSE(thetaMLE,OmegaSW,n)
	
	Results=list(inf=SEinf,sand=SEsw)
   	return(Results)
   	
}