/
betaNB-p-cor-nb.R
136 lines (136 loc) · 3.3 KB
/
betaNB-p-cor-nb.R
1
2
3
4
5
6
7
8
9
10
11
12
13
14
15
16
17
18
19
20
21
22
23
24
25
26
27
28
29
30
31
32
33
34
35
36
37
38
39
40
41
42
43
44
45
46
47
48
49
50
51
52
53
54
55
56
57
58
59
60
61
62
63
64
65
66
67
68
69
70
71
72
73
74
75
76
77
78
79
80
81
82
83
84
85
86
87
88
89
90
91
92
93
94
95
96
97
98
99
100
101
102
103
104
105
106
107
108
109
110
111
112
113
114
115
116
117
118
119
120
121
122
123
124
125
126
127
128
129
130
131
132
133
134
135
136
#' Estimate Squared Partial Correlation Coefficients
#' and Generate the Corresponding Sampling Distribution
#' Using Nonparametric Bootstrapping
#'
#' @author Ivan Jacob Agaloos Pesigan
#'
#' @details The vector of squared partial correlation coefficients
#' (\eqn{r^{2}_{p}})
#' is estimated from bootstrap samples.
#' Confidence intervals are generated by obtaining
#' percentiles corresponding to \eqn{100(1 - \alpha)\%}
#' from the generated sampling
#' distribution of \eqn{r^{2}_{p}},
#' where \eqn{\alpha} is the significance level.
#'
#' @return Returns an object
#' of class `betanb` which is a list with the following elements:
#' \describe{
#' \item{call}{Function call.}
#' \item{args}{Function arguments.}
#' \item{thetahatstar}{Sampling distribution of
#' \eqn{r^{2}_{p}}.}
#' \item{vcov}{Sampling variance-covariance matrix of
#' \eqn{r^{2}_{p}}.}
#' \item{est}{Vector of estimated
#' \eqn{r^{2}_{p}}.}
#' \item{fun}{Function used ("PCorNB").}
#' }
#'
#' @inheritParams BetaNB
#'
#' @examples
#' # Data ---------------------------------------------------------------------
#' data("nas1982", package = "betaNB")
#'
#' # Fit Model in lm ----------------------------------------------------------
#' object <- lm(QUALITY ~ NARTIC + PCTGRT + PCTSUPP, data = nas1982)
#'
#' # NB -----------------------------------------------------------------------
#' nb <- NB(
#' object,
#' R = 100, # use a large value e.g., 5000L for actual research
#' seed = 0508
#' )
#'
#' # PCorNB -------------------------------------------------------------------
#' out <- PCorNB(nb, alpha = 0.05)
#'
#' ## Methods -----------------------------------------------------------------
#' print(out)
#' summary(out)
#' coef(out)
#' vcov(out)
#' confint(out, level = 0.95)
#'
#' @family Beta Nonparametric Bootstrap Functions
#' @keywords betaNB pcor
#' @export
PCorNB <- function(object,
alpha = c(0.05, 0.01, 0.001)) {
stopifnot(
inherits(
object,
"nb"
)
)
if (object$lm_process$p < 2) {
stop("Two or more regressors is required.")
}
fun <- "PCorNB"
est <- .PCorSq(
srsq = .SPCor(
betastar = object$lm_process$betastar,
sigmacapx = object$lm_process$sigmacapx
)^2,
rsq = object$lm_process$rsq[1]
)
names(est) <- object$lm_process$xnames
foo <- function(x) {
sr <- .SPCor(
betastar = .BetaStarofSigma(
sigmacap = x,
q = 1 / sqrt(diag(x)),
k = object$lm_process$k
),
sigmacapx = x[
2:object$lm_process$k,
2:object$lm_process$k,
drop = FALSE
]
)
rsq <- .RSqofSigma(
sigmacap = x,
k = object$lm_process$k
)
return(
.PCorSq(
srsq = sr^2,
rsq = rsq
)
)
}
thetahatstar <- lapply(
X = object$thetahatstar,
FUN = foo
)
vcov <- stats::var(
do.call(
what = "rbind",
args = thetahatstar
)
)
colnames(vcov) <- rownames(vcov) <- names(est)
out <- list(
call = match.call(),
args = list(
object = object,
alpha = alpha
),
thetahatstar = thetahatstar,
jackknife = lapply(
X = object$jackknife,
FUN = foo
),
vcov = vcov,
est = est,
fun = fun
)
class(out) <- c(
"betanb",
class(out)
)
return(
out
)
}