[Distr-commits] r1274 - in branches/distr-2.8/pkg/distrMod: R man
noreply at r-forge.r-project.org
noreply at r-forge.r-project.org
Sat Aug 11 20:57:41 CEST 2018
Author: ruckdeschel
Date: 2018-08-11 20:57:40 +0200 (Sat, 11 Aug 2018)
New Revision: 1274
Modified:
branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R
branches/distr-2.8/pkg/distrMod/man/L2ParamFamily-class.Rd
branches/distr-2.8/pkg/distrMod/man/L2ParamFamily.Rd
Log:
[distrMod] branch 2.8:
+ more precise / explicite description of the requirements of slots L2deriv and L2deriv.fct in
the help files to generator L2ParamFamily and to L2ParamFamily-class.
+ the revised .CvMMDCovariance() uses vectorization in evaluation of random variables
and, wherever possible in integration; for the latter, this can be suppressed by
an argument useApply=TRUE through the ... argument
Modified: branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R
===================================================================
--- branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R 2018-08-11 14:33:45 UTC (rev 1273)
+++ branches/distr-2.8/pkg/distrMod/R/asCvMVarianceQtl.R 2018-08-11 18:57:40 UTC (rev 1274)
@@ -33,6 +33,8 @@
dotsInt[["upper"]] <- NULL
dotsInt[["stop.on.error"]] <- NULL
dotsInt[["distr"]] <- NULL
+ .useApply <- FALSE
+ if(!is.null(dotsInt$useApply)) .useApply <- dotsInt$useApply
if(missing(TruncQuantile)||TruncQuantile>1e-7) TruncQuantile <- 1e-8
@@ -147,11 +149,11 @@
if(is(distr,"AbscontDistribution")){
fqx <- function(x){qx <- q.l(distr)(x)
- return(sapply(qx,function(y)evalRandVar(L2deriv.0, y)))
+ evalRandVar(L2deriv.0, as.matrix(qx))[,,1]
}
- Delta0x.1 <- sapply(x.seq.1,fqx)
- Delta0x.2 <- sapply(x.seq.2,fqx)
- Delta0x.3 <- sapply(x.seq.3,fqx)
+ Delta0x.1 <- fqx(x.seq.1)
+ Delta0x.2 <- fqx(x.seq.2)
+ Delta0x.3 <- fqx(x.seq.3)
Delta0.1 <- h0/100*.csimpsum(Delta0x.1)
Delta0.2 <- rev(Delta0.1)[1]+h0*.csimpsum(Delta0x.2)
Delta0.3 <- rev(Delta0.2)[1]+h0/100*.csimpsum(Delta0x.3)
@@ -162,7 +164,7 @@
J <- do.call(myint, c(list(f=function(x) (Delta1.q(x)-J1)^2),dotsInt))
}else{
if(is(distr,"DiscreteDistribution")){
- L2x <- sapply(x.seq, function(x) evalRandVar(L2deriv.0, x))
+ L2x <- evalRandVar(L2deriv.0, as.matrix(x.seq))[,,1]
L2xdx <- L2x*prob
Delta0 <- cumsum(L2xdx)
J1 <- sum(Delta0*prob)
@@ -171,14 +173,14 @@
Delta.0 <- approxfun(x.seq, Delta, yleft = 0, yright = 0)
Delta <- Delta/J
}else{
- L2x <- function(x,y) (x<=y)*evalRandVar(L2deriv.0, x)
+ L2x <- function(x,y) (x<=y)*evalRandVar(L2deriv.0, as.matrix(x))[,,1]
Delta0 <- sapply(x.seq, function(Y){ fct <- function(x) L2x(x,y=Y)
- return(E(object=distr, fun = fct))})
+ return(E(object=distr, fun = fct, useApply = .useApply))})
Delta1 <- approxfun(x.seq, Delta0, yleft = 0, yright = 0)
Delta <- Delta1
- J1 <- E(object=distr, fun = Delta)
+ J1 <- E(object=distr, fun = Delta, useApply = .useApply)
Delta.0 <- function(x) Delta(x) - J1
- J <- E(object=distr, fun = function(x) Delta.0(x)^2 )
+ J <- E(object=distr, fun = function(x) Delta.0(x)^2, useApply = .useApply )
}
}
@@ -195,9 +197,9 @@
phiqx <- function(x){qx <- q.l(mu)(x)
return(phi(qx))}
- psi0qx.1 <- sapply(rev(x.mu.seq.1), phiqx)
- psi0qx.2 <- sapply(rev(x.mu.seq.2), phiqx)
- psi0qx.3 <- sapply(rev(x.mu.seq.3), phiqx)
+ psi0qx.1 <- phiqx(x.mu.seq.1)
+ psi0qx.2 <- phiqx(x.mu.seq.2)
+ psi0qx.3 <- phiqx(x.mu.seq.3)
psi0q.3 <- h0.mu/100*rev(.csimpsum(psi0qx.3))
psi0q.2 <- psi0q.3[1]+h0.mu*rev(.csimpsum(psi0qx.2))
psi0q.1 <- psi0q.2[1]+h0.mu/100*rev(.csimpsum(psi0qx.1))
@@ -207,13 +209,11 @@
psi.fct <- function(x) psi.q1(p(mu)(x))-psi1
}else{
if(is(mu,"DiscreteDistribution")&&is(distr,"DiscreteDistribution")){
+ Delta.mu <- phi(x.mu.seq)
+ pprob.mu <- p(distr)(x.mu.seq)
if(!all(support(mu)==support(distr))){
- Delta.mu <- sapply(x.mu.seq, phi)
- pprob.mu <- sapply(x.mu.seq, p(distr))
- L2x.mu <- sapply(x.mu.seq, function(x) evalRandVar(L2deriv.0, x))
+ L2x.mu <- evalRandVar(L2deriv.0, as.matrix(x.mu.seq))[,,1]
}else{
- Delta.mu <- sapply(x.mu.seq, phi)
- pprob.mu <- cumsum(prob)
L2x.mu <- L2x
}
psi1 <- sum(pprob.mu*Delta.mu*prob.mu)
@@ -223,11 +223,11 @@
}else{
## integrand phi x Ptheta in formula (51) [ibid]
phi1 <- function(x) phi(x) * p(distr)(x)
- psi1 <- E(object = mu, fun = phi1)
+ psi1 <- E(object = mu, fun = phi1, useApply = .useApply)
phixy <- function(x,y) (x<=y)*phi(y)
psi0 <- sapply(x.mu.seq, function(X){ fct <- function(y) phixy(x=X,y=y)
- return(E(object=mu, fun = fct))})
+ return(E(object=mu, fun = fct, useApply = .useApply))})
psi.1 <- approxfun(x.mu.seq, psi0, yleft = 0, yright = rev(psi0)[1])
psi.fct <- function(x) psi.1(x)-psi1
}
@@ -240,8 +240,7 @@
E1 <- do.call(myint, c(list(f=psi.q),dotsInt))
E3 <- do.call(myint, c(list(f=function(x){
qx <- q.l(distr)(x)
- L2qx <- sapply(qx,function(y)
- evalRandVar(L2deriv.0, y))
+ L2qx <- evalRandVar(L2deriv.0,as.matrix(qx))[,,1]
return(psi.fct(qx)*L2qx)
}), dotsInt))
psi.01.f <- function(x) (psi.fct(x)-E1)/E3
@@ -250,8 +249,7 @@
if(is(distr,"DiscreteDistribution")){
## E2 = Cov_mu (psi)
# E2 <- sum(psi0^2*prob)
- psi0 <- sapply(x.seq, psi.fct)
-
+ psi0 <- psi.fct(x.seq)
E1 <- sum(psi0*prob)
E3 <- sum(psi0*L2x*prob)
psi.01d <- (psi0-E1)/E3
@@ -259,12 +257,13 @@
psi.01.f <- function(x) (psi.fct(x)-E1)/E3*liesInSupport(distr,x)
}else{
## E2 = Cov_mu (psi)
-# E2 <- E(object=distr, fun = function(x) psi(x)^2)
- L2x <- function(x,y) (x<=y)*evalRandVar(L2deriv.0, x)
- E1 <- E(object=distr, fun = psi.fct )
- E3 <- E(object=distr, fun = function(x) psi.fct(x)*evalRandVar(L2deriv.0, x))
+# E2 <- E(object=distr, fun = function(x) psi(x)^2, useApply = .useApply)
+ L2x <- function(x,y) (x<=y)*evalRandVar(L2deriv.0, as.matrix(x))[,,1]
+ E1 <- E(object=distr, fun = psi.fct, useApply = .useApply )
+ E3 <- E(object=distr, fun = function(x)
+ psi.fct(x)*evalRandVar(L2deriv.0, as.matrix(x))[,,1], useApply = .useApply)
psi.01.f <- function(x) (psi.fct(x) - E1)/E3
- E4 <- E(object=distr, fun = function(x) psi.01.f(x)^2)
+ E4 <- E(object=distr, fun = function(x) psi.01.f(x)^2, useApply = .useApply)
}
}
@@ -293,9 +292,10 @@
for(i in 1:Dim)
{ if(is(distr,"AbscontDistribution")){
- fct0.q1 <- sapply(x.seq.1, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
- fct0.q2 <- sapply(x.seq.2, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
- fct0.q3 <- sapply(x.seq.3, function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))})
+ fct.q <- function(x){qx <- q.l(distr)(x); return(L2deriv.0 at Map[[i]](qx))}
+ fct0.q1 <- fct.q(x.seq.1)
+ fct0.q2 <- fct.q(x.seq.2)
+ fct0.q3 <- fct.q(x.seq.3)
Delta0.q1 <- h0/100*.csimpsum(fct0.q1)
Delta0.q2 <- rev(Delta0.q1)[1]+h0*.csimpsum(fct0.q2)
Delta0.q3 <- rev(Delta0.q2)[1]+h0/100*.csimpsum(fct0.q3)
@@ -307,7 +307,7 @@
assign("Delta1.q", Delta1.q, envir=env.i)
}else{
if(is(distr,"DiscreteDistribution")){
- L2x <- sapply(x.seq, function(x) evalRandVar(L2deriv.0, x)[i])
+ L2x <- evalRandVar(L2deriv.0, as.matrix(x.seq))[i,,1]
L2xdx <- L2x*prob
Delta.0 <- cumsum(L2xdx)
Delta.f <- approxfun(x.seq, Delta.0, yleft = 0, yright = 0)
@@ -318,7 +318,7 @@
}else{
fct0 <- function(x,y) L2deriv.0 at Map[[i]](x)*(x<=y)
Delta0 <- sapply(x.seq, function(Y){ fct <- function(x) fct0(x,y=Y)
- return(E(object=distr, fun = fct))})
+ return(E(object=distr, fun = fct, useApply=.useApply))})
Delta1 <- approxfun(x.seq, Delta0, yleft = 0, yright = 0)
if(is(distr,"DiscreteDistribution"))
Delta <- function(x) Delta1(x) * (x %in% support(distr))
@@ -350,12 +350,12 @@
## J = Var_Ptheta Delta
##-t-## print(system.time({
- J1 <- E(object=distr, fun = Delta)
+ J1 <- E(object=distr, fun = Delta)#, useApply = .useApply)
##-t-## }))
Delta.0 <- Delta - J1
##-t-## print(system.time({
- J <- E(object=distr, fun = Delta.0 %*%t(Delta.0))
+ J <- E(object=distr, fun = Delta.0 %*%t(Delta.0))#, useApply = .useApply)
##-t-## }))
### CvM-IC phi
phi <- as(distr::solve(J)%*%Delta.0,"EuclRandVariable")
@@ -371,7 +371,7 @@
phi1 <- EuclRandVariable(Map = Map.phi1, Domain = Reals())
##-t-## print(system.time({
- psi1 <- E(object=mu, fun = phi1)
+ psi1 <- E(object=mu, fun = phi1)#, useApply = .useApply)
##-t-## }))
## obtaining IC psi (formula (51))
@@ -389,7 +389,7 @@
qxm <- q.l(mu)(x.mu.seq.b)
##-t-## print(system.time({
- fct0.qq <- sapply(qxm, phi at Map[[i]])
+ fct0.qq <- phi at Map[[i]](qxm)
##-t-## }))
fct0.q1 <- rev(fct0.qq[iN.mu.1])
fct0.q2 <- rev(fct0.qq[iN.mu.2])
@@ -405,7 +405,7 @@
assign("psi0", psi0, envir=env.i)
}else{
if(is(mu,"DiscreteDistribution")){
- phi.mu <- sapply(x.mu.seq, function(x) evalRandVar(phi,x)[i])
+ phi.mu <- evalRandVar(phi,as.matrix(x.mu.seq))[i,,1]
psi0.d <- cumsum(phi.mu*prob.mu)
psi0.a <- approxfun(x.mu.seq, psi0.d, yleft = 0, yright = 0)
psi0 <- function(x) psi0.a(x)
@@ -414,11 +414,11 @@
assign("psi0.d", psi0.d, envir=env.i)
assign("psi0", psi0, envir=env.i)
}else{
- fct0 <- function(x,y) evalRandVar(phi, y)[i]*(x<=y)
+ fct0 <- function(x,y) evalRandVar(phi, as.matrix(y))[i,,1]*(x<=y)
phi0 <- sapply(x.mu.seq,
function(X){
fct <- function(y) fct0(x = X, y)
- return(E(object = mu, fun = fct))
+ return(E(object = mu, fun = fct, useApply = .useApply))
})
phi0a <- approxfun(x.mu.seq, phi0, yleft = 0, yright = rev(phi0)[1])
if(is(distr,"DiscreteDistribution"))
@@ -441,13 +441,13 @@
# print(Map.psi)
psi <- EuclRandVariable(Map = Map.psi, Domain = Reals())
-# E2 <- E(object=distr, fun = psi %*%t(psi))
+# E2 <- E(object=distr, fun = psi %*%t(psi), useApply = .useApply)
## E2 = Cov_mu (psi)
### control: centering & standardization
L2deriv.0 <- L2Fam at L2deriv[[1]]
##-t-## print(system.time({
- E1 <- E(object=distr, fun = psi )
+ E1 <- E(object=distr, fun = psi)
##-t-## }))
##-t-## print(system.time({
E3 <- E(object=distr, fun = psi %*%t(L2deriv.0))
@@ -563,8 +563,6 @@
}
L2deriv <- L2deriv(L2Fam)[[1]]
-# y.seq <- sapply(x.seq, function(x) evalRandVar(L2deriv, x))
-# plot(x.seq[!is.na(y.seq)],y.seq ,type="l")
## are we working with a one-dim L2deriv or not?
Modified: branches/distr-2.8/pkg/distrMod/man/L2ParamFamily-class.Rd
===================================================================
--- branches/distr-2.8/pkg/distrMod/man/L2ParamFamily-class.Rd 2018-08-11 14:33:45 UTC (rev 1273)
+++ branches/distr-2.8/pkg/distrMod/man/L2ParamFamily-class.Rd 2018-08-11 18:57:40 UTC (rev 1274)
@@ -69,7 +69,10 @@
properties of the family. }
\item{\code{L2deriv}}{
object of class \code{"EuclRandVariable"}:
- L2 derivative of the family. }
+ L2 derivative of the family. Its map slot must contain a list of functions.
+ Each function in this list must have just one argument \code{x},
+ which is vectorized, (i.e., callable for a vector-valued \code{x}),
+ and has a one-dimensional, numeric return value.}
\item{\code{L2deriv.fct}}{
object of class \code{"function"}: mapping from
the parameter space (argument \code{param} of class
@@ -77,18 +80,32 @@
value of the L2derivative; \code{L2deriv.fct} is then used from observation
\code{x} to value of the L2derivative; \code{L2deriv.fct} is used by
\code{modifyModel} to move the L2deriv according to a change in the
- parameter }
- \item{\code{L2derivSymm}}{[
- object of class \code{"FunSymmList"}:
- symmetry of the maps included in \code{L2deriv}. }
- \item{\code{L2derivDistr}}{
- object of class \code{"OptionalDistrListOrCall"} (i.e., \code{NULL} or
- an object of class \code{"DistrList"} or the respective call to generate
- the latter object): if non-null and non-call, a
- list which includes the distribution of \code{L2deriv}. }
- \item{\code{L2derivDistrSymm}}{
- object of class \code{"DistrSymmList"}:
- symmetry of the distributions included in \code{L2derivDistr}. }
+ parameter.
+ More specifically, let us call the parts \code{main} and \code{nuisance}
+ of the parameter the \emph{unknown} parameter. If this unknown parameter is
+ one-dimensional, the return value of \code{L2deriv.fct} must be a function
+ in argument \code{x}, which is vectorized, (i.e.,
+ callable for a vector-valued \code{x}), and has a one-dimensional, numeric
+ return value. In case the dimension of the unknown parameter is larger
+ than one, the return value must be a list of functions, each of which
+ satisfies the conditions formulated for the case of a one-dimensional
+ parameter of interest. The order of the components of this list is
+ the same as the order of the parameter coordinates in \code{main}, followed
+ by the ones in \code{nuisance}.}
+ \item{L2derivSymm}{ object of class \code{"FunSymmList"}:
+ symmetry of the maps contained in \code{L2deriv}; a list
+ of symmetry properties of the same length as the return value of
+ \code{ L2deriv.fct }.}
+ \item{L2derivDistr}{object of class \code{"OptionalDistrListOrCall"}
+ (i.e., \code{NULL} or an object of class \code{"DistrList"} or
+ the respective call to generate the latter object): if non-null
+ and non-call, a list which includes the distribution of \code{L2deriv};
+ the length of this list of univariate distributions must be of the same
+ length as the return value of \code{ L2deriv.fct }.}
+ \item{L2derivDistrSymm}{ object of class \code{"DistrSymmList"}:
+ symmetry of the distributions contained in \code{L2derivDistr};
+ the length of this list of symmetry properties must be
+ of the same length as the return value of \code{ L2deriv.fct }. }
\item{\code{FisherInfo.fct}}{
object of class \code{"function"}:
mapping from the parameter space (argument \code{param} of class
Modified: branches/distr-2.8/pkg/distrMod/man/L2ParamFamily.Rd
===================================================================
--- branches/distr-2.8/pkg/distrMod/man/L2ParamFamily.Rd 2018-08-11 14:33:45 UTC (rev 1273)
+++ branches/distr-2.8/pkg/distrMod/man/L2ParamFamily.Rd 2018-08-11 18:57:40 UTC (rev 1274)
@@ -51,13 +51,31 @@
\code{param} of class \code{"ParamFamParameter"}) to a mapping from
observation \code{x} to the value of the L2derivative;
\code{L2deriv.fct} is used by \code{modifyModel} to
- move the L2deriv according to a change in the parameter }
+ move the L2deriv according to a change in the parameter,
+ and to fill slot \code{L2deriv}.
+ More specifically, let us call the parts \code{main} and \code{nuisance}
+ of the parameter the \emph{unknown} parameter. If this unknown parameter is
+ one-dimensional, the return value of \code{L2deriv.fct} must be a function
+ in argument \code{x}, which is vectorized, (i.e.,
+ callable for a vector-valued \code{x}), and has a one-dimensional, numeric
+ return value. In case the dimension of the unknown parameter is larger
+ than one, the return value must be a list of functions, each of which
+ satisfies the conditions formulated for the case of a one-dimensional
+ parameter of interest. The order of the components of this list is
+ the same as the order of the parameter coordinates in \code{main}, followed
+ by the ones in \code{nuisance}.}
\item{L2derivSymm}{ object of class \code{"FunSymmList"}:
- symmetry of the maps contained in \code{L2deriv} }
+ symmetry of the maps contained in \code{L2deriv}; a list
+ of symmetry properties of the same length as the return value of
+ \code{ L2deriv.fct }.}
\item{L2derivDistr}{ object of class \code{"UnivarDistrList"}:
- distribution of \code{L2deriv} }
+ distribution of \code{L2deriv}; the length of this list
+ of univariate distributions must be of the same length as the
+ return value of \code{ L2deriv.fct }.}
\item{L2derivDistrSymm}{ object of class \code{"DistrSymmList"}:
- symmetry of the distributions contained in \code{L2derivDistr} }
+ symmetry of the distributions contained in \code{L2derivDistr};
+ the length of this list of symmetry properties must be
+ of the same length as the return value of \code{ L2deriv.fct }. }
\item{FisherInfo.fct}{function: mapping from the parameter space (argument
\code{param} of class \code{"ParamFamParameter"}) to the set of positive
semidefinite matrices; \code{FisherInfo.fct} is used by \code{modifyModel} to
More information about the Distr-commits
mailing list