copula (version 0.999-19.1)

pobs: Pseudo-Observations

Description

Compute the pseudo-observations for the given data matrix.

Usage

pobs(x, na.last = "keep",
     ties.method = eval(formals(rank)$ties.method), lower.tail = TRUE)

Arguments

x

\(n\times d\)-matrix (or \(d\)-vector) of random variates to be converted to pseudo-observations.

na.last

string passed to rank; see there.

ties.method

string specifying how ranks should be computed if there are ties in any of the coordinate samples of x; passed to rank.

lower.tail

logical which, if FALSE, returns the pseudo-observations when applying the empirical marginal survival functions.

Value

matrix (or vector) of the same dimensions as x containing the pseudo-observations.

Details

Given \(n\) realizations \(\bm{x}_i=(x_{i1},\dots,x_{id})^T\), \(i\in\{1,\dots,n\}\) of a random vector \(\bm{X}\), the pseudo-observations are defined via \(u_{ij}=r_{ij}/(n+1)\) for \(i\in\{1,\dots,n\}\) and \(j\in\{1,\dots,d\}\), where \(r_{ij}\) denotes the rank of \(x_{ij}\) among all \(x_{kj}\), \(k\in\{1,\dots,n\}\). When there are no ties in any of the coordinate samples of x, the pseudo-observations can thus also be computed by component-wise applying the marginal empirical distribution functions to the data and scaling the result by \(n/(n+1)\). This asymptotically negligible scaling factor is used to force the variates to fall inside the open unit hypercube, for example, to avoid problems with density evaluation at the boundaries. Note that pobs(, lower.tail=FALSE) simply returns 1-pobs().

Examples

# NOT RUN {
## Simple definition of the function:
pobs
# }
# NOT RUN {
<!-- %packageDescription("Matrix")% for debugging checks / Matrix -->
# }
# NOT RUN {
## Draw from a multivariate normal distribution
d <- 10
set.seed(1)
P <- Matrix::nearPD(matrix(pmin(pmax(runif(d*d), 0.3), 0.99), ncol=d))$mat
diag(P) <- rep(1, d)
n <- 500
x <- MASS::mvrnorm(n, mu = rep(0, d), Sigma = P)

## Compute pseudo-observations (should roughly follow a Gauss
## copula with correlation matrix P)
u <- pobs(x)
plot(u[,5],u[,10], xlab=quote(italic(U)[1]), ylab=quote(italic(U)[2]))
# }
# NOT RUN {
## All components: pairwise plot
pairs(u, gap=0, pch=".", labels =
      as.expression( lapply(1:d, function(j) bquote(italic(U[.(j)]))) ))
# }