Nothing
TMLG <-
function(X,y,delta, Beta.t,sigma.t, const,tu,control) {
# TML, Gaussian case, iteratively reweighting algorithm
maxit.sigma <- control$maxit.sigma
tol.sigma <- control$tol.sigma
maxit.Beta <- control$maxit.Beta
tol.Beta <- control$tol.Beta
Maxit.TML <- control$Maxit.TML
tol.TML.sigma <- control$tol.TML.sigma
tol.TML.Beta <- control$tol.TML.Beta
alg.sigma <- control$alg.sigma
nitmon <- control$nitmon
p <- ncol(X); Nit <- 1; Mxt <- 10
Beta0 <- Beta.t
sigma0 <- sigma.t
mui.t <- X %*% as.matrix(Beta.t)
rs0 <- (y-mui.t)/sigma.t
wgt <- ww(rs0,-tu, tu)
DeltaB <- Deltas <- 100000
Nit <- 1
while ( ( Deltas > tol.TML.sigma | DeltaB > tol.TML.Beta ) & Nit <= Maxit.TML) {
mui <- X %*% as.matrix(Beta0) # check
if (alg.sigma ==1) {
zs <- TML.SigmaG(X,y,delta, sigma0, sigma.t, mui,mui.t, wgt,const,tu, tol=tol.sigma,maxit=maxit.sigma,nitmon)
sigma <- zs$sigma; nits <- zs$nit }
if (alg.sigma !=1) {
# determine lower and upper values for sigma
smu <- sml <- sigma0; nit <- 1
funu <- funl <- TML.Ave2G(X,y,delta, smu,sigma.t, mui,mui.t, wgt,tu)-const
while(funu > 0 & nit <= Mxt) {
smu <- smu*1.5
funu <- TML.Ave2G(X,y,delta, smu,sigma.t, mui,mui.t, wgt,tu)-const
nit <- nit+1}
if (nitmon) cat("smu, funu, nit",smu,funu,nit,"\n")
nit <- 1
while(funl < 0 & nit <= Mxt) {
sml <- sml*0.5
funl <- TML.Ave2G(X,y,delta, sml,sigma.t, mui,mui.t, wgt,tu)-const
nit <- nit+1}
if (nitmon) cat("sml, funl, nit",sml,funl,nit,"\n")
# compute new sigma
zs <- regfal(TML.Ave2G,cc=const,lower=sml,upper=smu,nint=5,tol=tol.sigma,maxit=maxit.sigma,
X=X,y=y,delta=delta, sigma.t=sigma.t, mui=mui, mui.t=mui.t, wgt=wgt,cu=tu)
sigma <- zs$solution; nits <- zs$nit}
if (nitmon) cat("sigma", nits, "\n")
# compute new Beta
zB <- TML.BetaG(X,y,delta,Beta0,sigma, Beta.t,sigma.t, tu, maxit=maxit.Beta,tol=tol.Beta,nitmon)
Beta <- zB$Beta
DeltaB <- max(abs(Beta - Beta0)); Deltas <- abs(sigma-sigma0)
if (nitmon) cat("Nit=",Nit,"DeltaB",DeltaB,"Deltas",Deltas, "\n")
Beta0 <- Beta; sigma0 <- sigma
Nit <- Nit+1}
zres <- list(Beta=Beta,sigma=sigma,Nit=Nit)
zres}
Any scripts or data that you put into this service are public.
Add the following code to your website.
For more information on customizing the embed code, read Embedding Snippets.