A vectorised uniroot function would be useful for function inversion,
e.g. for quantile functions and random number generation. To address
this, I have implemented rstpm2::vuniroot that adapts the C function
R_zeroin2 for Brent's method for a vectorised objective. The function
currently uses Rcpp, but could be re-implemented using the C API.
As an example, we can now rapidly sample from a proportional hazards
mixture Weibull distribution:
pweibullMixturePH <- function(q, p1, RR, shape1, shape2, scale1=1, scale2=1)
1 - (p1*pweibull(q, shape1, scale1, lower.tail=FALSE) +
(1-p1)*pweibull(q, shape2, scale2, lower.tail=FALSE))^RR
rfun <- function(pfun) function(n, ..., lower=1e-5, upper=1e6) {
u <- runif(n)
objective <- function(q) pfun(q, ...) - u
rstpm2::vuniroot(objective, lower=rep(lower,length=n),
upper=rep(upper,length=n))$root
}
rweibullMixturePH <- rfun(pweibullMixturePH)
set.seed(12345)
y <- rweibullMixturePH(n=1e4,p1=0.5,RR=2,shape1=1.5,shape2=0.5)
Has anyone previously developed a similar vectorised uniroot function?
Finally, would this be a useful addition to core R?
-- Mark
När du skickar e-post till Karolinska Institutet (KI) innebär detta att KI
kommer att behandla dina personuppgifter. Här finns information om hur KI
behandlar personuppgifter<https://ki.se/medarbetare/integritetsskyddspolicy>.
Sending email to Karolinska Institutet (KI) will result in KI processing your
personal data. You can read more about KI’s processing of personal data
here<https://ki.se/en/staff/data-protection-policy>.
______________________________________________
R-help@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.
