# R/utils.R In popKorn: For interval estimation of mean of selected populations

#### Defines functions computeMaxintegrate4bintegrate4aintegrate2bintegrate2aintegrand4bintegrand4aintegrand2bintegrand2acomponent3component2component1bcomponent1ah2optimalCScalar2

```# Utility functions that are not exported
optimalCScalar2 <- function(lambda, alpha, min.loc, n, p, k, var.known) {
out <- uniroot(f=eqnForMin, interval=c(0,500), lambda=lambda,
alpha=alpha, min.loc=min.loc, n=n, p=p, k=k, var.known=var.known)\$root

return(out)
}

# For finding optimal lambda
h2 <- function(lambda, alpha, n, p, k, var.known) {
abs(optimalC(lambda, alpha, 'zero', n, p, k, var.known) -
optimalC(lambda, alpha, 'infty', n, p, k, var.known))
}
h2.v <- Vectorize(h2, vectorize.args="lambda")

####### Unknown variance case, k=1 #########################
# For the unknown variance, single p case
component1a <- function(x, c.val, lambda, p)
pnorm(c.val*x)^p - pnorm(-lambda*c.val*x)^p
component1b <- function(x, c.val, lambda)
pnorm(c.val*x) - pnorm(-lambda*c.val*x)

component2 <- function(x, c.val, lambda, n, p) {
dgamma(x, shape=0.5*p*(n-1), scale=2/(p*(n-1)))
}

component3 <- function(x, c.val, lambda, n, p) {
log.scale <- 0.5*p*(n-1)*log(x) + 0.5*x*p*(n-1) - 0.5*p*(n-1)*x^2
2*exp(log.scale)
}

#2a corresponds to 'zero'
#2b corresponds to 'infty'
integrand2a <- function(x, c.val, lambda, n, p) {
component1a(x,c.val,lambda,p)*
component2(x,c.val,lambda,n,p)*
component3(x,c.val,lambda,n,p)
}
integrand2b <- function(x, c.val, lambda, n, p) {
component1b(x,c.val,lambda)*
component2(x,c.val,lambda,n,p)*
component3(x,c.val,lambda,n,p)
}

####### Unknown variance case, k>1 #########################
# For the unknown variance, multiple k case

#4a corresponds to 'zero'
#4b corresponds to 'infty'
integrand4a <- function(x, c.val, lambda, n, p, k) {
component1a(x,c.val,lambda,p-k+1)*
component2(x,c.val,lambda,n,p)*
component3(x,c.val,lambda,n,p)*
component1b(x,c.val,lambda)^(k-1)
}
integrand4b <- function(x, c.val, lambda, n, p, k) {
(component1b(x,c.val,lambda)^k)*
component2(x,c.val,lambda,n,p)*
component3(x,c.val,lambda,n,p)
}

###### integration functions ###########################

## Case k=1, min. at 0
integrate2a <- function(lambda, c.val, n, p) {
integrate(integrand2a, 0, Inf, c.val=c.val, lambda=lambda, n=n, p=p)\$value
}

## Case k=1, min. at infty
integrate2b <- function(lambda, c.val, n, p) {
integrate(integrand2b, 0, Inf, c.val=c.val, lambda=lambda, n=n, p=p)\$value
}

## Case k>1, min. at zero
integrate4a <- function(lambda, c.val, n, p, k) {
integrate(integrand4a, 0, Inf, c.val=c.val, lambda=lambda, n=n, p=p, k=k)\$value
}

## Case k>1, min. at infty
integrate4b <- function(lambda, c.val, n, p, k) {
integrate(integrand4b, 0, Inf, c.val=c.val, lambda=lambda, n=n, p=p, k=k)\$value
}

####### statistic() function for bootstrapping ############
computeMax <- function(d, i, k=1) {
d.2 <- d\$values[i]
sample.means <- tapply(d.2, d\$ind, mean)
sample.means <- sort(sample.means, decreasing=TRUE)
sample.means[1:k]
}
```

## Try the popKorn package in your browser

Any scripts or data that you put into this service are public.

popKorn documentation built on May 29, 2017, 10:32 p.m.