Nothing
#if(FALSE){
## faster check for ContICs
setMethod("checkIC", signature(IC = "ContIC", L2Fam = "L2ParamFamily"),
function(IC, L2Fam, out = TRUE, forceContICMethod = FALSE, ..., diagnostic = FALSE){
D1 <- L2Fam@distribution
if( dimension(Domain(IC@Curve[[1]])) != dimension(img(D1)))
stop("dimension of 'Domain' of 'Curve' != dimension of 'img' of 'distribution' of 'L2Fam'")
res <- .prepareCheckMakeIC(L2Fam, w = IC@weight, forceContICMethod, ..., diagnostic = diagnostic)
## if it pays off to use symmetry/ to compute integrals in L2deriv space
## we compute the following integrals:
## G1 = E w, G2 = E Lambda w, G3 = E Lambda Lambda' w
## we want to compute:
## Delta1 = E (A Lambda-a) w, Delta2 = E (A Lambda-a) Lambda' w
## where A = stand(IC), a=cent(IC)
## hence Delta1 = A G2 - a G1, Delta2 = A G3 - a G2'
### otherwise the return value is NULL and we use the standard method
if(is.null(res))
return(getMethod("checkIC", signature(IC = "IC",
L2Fam = "L2ParamFamily"))(IC,L2Fam, out = out, ..., diagnostic = diagnostic))
A <- stand(IC); a <- cent(IC)
G1 <- res$G1; G2 <- res$G2; G3 <- res$G3
Delta1 <- A%*%G2- a*G1
Delta2 <- A%*%G3 - a%*%t(G2)
trafoL <- trafo(L2Fam)
Delta2 <- Delta2 - trafoL
Prec <- ceiling(12-round(max(log(abs(trafoL)+1e-14,10)))/2)
## PR 20190407: in output in if(out)
## deleting all digits beyond 1e-12 (as numeric fuzz) --
## but check for relative accuracy by means of the "size" of the Fisher information
## measured in by the max(trafo)
if(out){
cent.out <- round(Delta1*10^Prec)/10^Prec
cat("precision of centering:\t", cent.out, "\n")
oldOps <- options()
on.exit(do.call(options,oldOps))
consist.out <- round(Delta2*10^Prec)/10^Prec
options(digits=5,scipen=-2)
cat("precision of Fisher information:\n")
print(consist.out)
cat("precision of Fisher information - relativ error [%]:\n")
relconsist.out <- round(Delta2/trafoL*10^(Prec+2))/10^Prec
class(relconsist.out) <- c("relMatrix",class(consist.out))
print(relconsist.out)
if(diagnostic){
print(attr(res$G1, "diagnostic"))
print(attr(res$G2, "diagnostic"))
print(attr(res$G3, "diagnostic"))
}
}
prec <- max(abs(Delta1), abs(Delta2))
names(prec) <- "maximum deviation"
if(diagnostic) attr(prec, "diagnostic") <- c(attr(res$G1, "diagnostic"),
attr(res$G2, "diagnostic"), attr(res$G3, "diagnostic"))
return(prec)
})
## make some L2function a pIC at a model
setMethod("makeIC", signature(IC = "ContIC", L2Fam = "L2ParamFamily"),
function(IC, L2Fam, forceContICMethod = FALSE, ..., diagnostic = FALSE){
D1 <- L2Fam@distribution
if( dimension(Domain(IC@Curve[[1]])) != dimension(img(D1)))
stop("dimension of 'Domain' of 'Curve' != dimension of 'img' of 'distribution' of 'L2Fam'")
dims <- nrow(trafo(L2Fam))
if(dimension(IC@Curve) != dims)
stop("Dimension of IC and parameter must be equal")
res <- .prepareCheckMakeIC(L2Fam, w = IC@weight, forceContICMethod, ..., diagnostic = diagnostic)
if(diagnostic &&!is.null(res)){
print(attr(res$G1, "diagnostic"))
print(attr(res$G2, "diagnostic"))
print(attr(res$G3, "diagnostic"))
}
## if it pays off to use symmetry/ to compute integrals in L2deriv space
## we compute the following integrals:
## G1 = E w, G2 = E Lambda w, G3 = E Lambda Lambda' w
## we want to compute:
## Delta1 = E (A Lambda-a) w, Delta2 = E (A Lambda-a) Lambda' w
## where A = stand(IC), a=cent(IC)
## hence Delta1 = A G2 - a G1, Delta2 = A G3 - a G2'
### otherwise the return value is NULL and we use the standard method
if(is.null(res))
return(getMethod("makeIC", signature(IC = "IC",
L2Fam = "L2ParamFamily"))(IC,L2Fam,..., diagnostic = diagnostic))
G1 <- res$G1; G2 <- res$G2; G3 <- res$G3
trafO <- trafo(L2Fam@param)
nrvalues <- nrow(trafO)
dims <- ncol(trafO)
cent0 <- c(G2/G1)
stand1 <- trafO%*%distr::solve(G3-cent0%*%t(G2))
cent1 <- c(stand1%*%cent0)
# print(list(stand1,stand(IC),cent1,cent(IC)))
L2.f <- as(diag(nrvalues) %*% L2Fam@L2deriv , "EuclRandVariable")
D1 <- L2Fam@distribution
IC1.f <- function(x){ indS <- liesInSupport(D1,x,checkFin=TRUE)
Lx <- sapply(x, function(y) evalRandVar(L2.f,y))
indS* (stand1%*%Lx-cent1) * weight(IC@weight)(Lx)}
IC1.l <- vector("list",nrvalues)
for(i in 1:nrvalues){
IC1.l[[i]] <- function(x){}
body(IC1.l[[i]]) <- substitute( c((IC1.s(x))[i,]), list(IC1.s=IC1.f, i=i))
}
IC1.c <- EuclRandVariable(Map = IC1.l, Domain = Domain(IC@Curve[[1]]),
Range = Reals())
cIC1 <- new("ContIC")
cIC1@name <- IC@name
cIC1@Curve <- EuclRandVarList(IC1.c)
cIC1@Risks <- IC@Risks
cIC1@Infos <- IC@Infos
cIC1@CallL2Fam <- L2Fam@fam.call
cIC1@clip <- IC@clip
cIC1@cent <- cent1
cIC1@stand <- stand1
cIC1@lowerCase <- IC@lowerCase
cIC1@neighborRadius <- IC@neighborRadius
cIC1@weight <- IC@weight
cIC1@biastype <- IC@biastype
cIC1@normtype <- IC@normtype
cIC1@modifyIC <- IC@modifyIC
addInfo(cIC1) <- c("IC<-",
"generated by affine linear trafo to enforce consistency")
if(diagnostic) attr(cIC1, "diagnostic") <- c(attr(res$G1, "diagnostic"),
attr(res$G2, "diagnostic"), attr(res$G3, "diagnostic"))
return(cIC1)
})
.prepareCheckMakeIC <- function(L2Fam, w, forceContICMethod, ..., diagnostic = FALSE){
dims <- length(L2Fam@param)
trafo <- trafo(L2Fam@param)
nrvalues <- nrow(trafo)
z.comp <- rep(TRUE,dims)
A.comp <- matrix(rep(TRUE,dims^2),nrow=dims)
to.comp.i <- (dims+1)*(dims+2)/2
to.comp.a <- (dims+1)*nrvalues
L2deriv <- as(diag(dims) %*% L2Fam@L2deriv, "EuclRandVariable")
z.comp <- rep(TRUE,dims)
A.comp <- matrix(TRUE, dims, dims)
# print(list(z.comp,A.comp))
# otherwise if nrvalues > 1 # formerly: trafo == unitMatrix #
# may use symmetry info
if(dims>1){
comp <- .getComp(L2deriv, L2Fam@distrSymm, L2Fam@L2derivSymm, L2Fam@L2derivDistrSymm)
z.comp <- comp$"z.comp"
A.comp <- comp$"A.comp"
t.comp.i <- sum(z.comp)+sum(A.comp)+1
}
if(to.comp.a < to.comp.i && !forceContICMethod) return(NULL)
res <- .getG1G2G3Stand(L2deriv = L2deriv, Distr = L2Fam@distribution,
A.comp = A.comp, z.comp = z.comp, w = w, ...,
diagnostic = diagnostic)
return(res)
}
.getG1G2G3Stand <- function(L2deriv, Distr, A.comp, z.comp, w, ..., diagnostic = FALSE){
dotsI <- .filterEargs(list(...))
if(is.null(dotsI$useApply)) dotsI$useApply <- FALSE
w.fct <- function(x){
weight(w)(evalRandVar(L2deriv, as.matrix(x)) [,,1])
}
integrand2 <- function(x, L2.i){
return(L2.i(x)*w.fct(x))
}
diagn <- if(diagnostic) vector("list", sum(z.comp)+sum(A.comp))
if(diagnostic) dotsI$diagnostic <- TRUE
Eargs <- c(list(object = Distr, fun = w.fct), dotsI)
res1 <- do.call(E,Eargs)
k <- 0
nrvalues <- length(L2deriv)
res2 <- numeric(nrvalues)
for(i in 1:nrvalues){
if(z.comp[i]){
integrand2i <- function(x) integrand2(x,L2deriv@Map[[i]])
Eargs <- c(list(object = Distr, fun = integrand2i), dotsI)
res2[i] <- buf <- do.call(E,Eargs)
if(diagnostic){k <- k + 1; diagn[[k]] <- attr(buf,"diagnostic")}
}else{
res2[i] <- 0
}
}
if(diagnostic) {k1 <- k; attr(res2, "diagnostic") <- diagn[(1:k1)]}
cent <- res2/res1
integrandA <- function(x, L2.i, L2.j, i, j){
return((L2.i(x) - cent[i])*(L2.j(x) - cent[j])*w.fct(x = x))
}
nrvalues <- length(L2deriv)
erg <- matrix(0, ncol = nrvalues, nrow = nrvalues)
for(i in 1:nrvalues){
for(j in i:nrvalues){
if(A.comp[i,j]){
integrandAij <- function(x) integrandA(x, L2.i = L2deriv@Map[[i]],
L2.j = L2deriv@Map[[j]], i = i, j = j)
Eargs <- c(list(object = Distr, fun = integrandAij), dotsI)
erg[i, j] <- buf <- do.call(E,Eargs)
if(diagnostic){k <- k + 1; diagn[[k]] <- attr(buf,"diagnostic")}
}
}
}
erg[col(erg) < row(erg)] <- t(erg)[col(erg) < row(erg)]
if(diagnostic) {k1 <- k; attr(erg, "diagnostic") <- diagn[-(1:k1)]}
return(list(G1=res1,G2=res2, G3=erg))
}
#}
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.