Last updated on 2019-12-05 00:51:48 CET.
| Flavor | Version | Tinstall | Tcheck | Ttotal | Status | Flags |
|---|---|---|---|---|---|---|
| r-devel-linux-x86_64-debian-clang | 0.3.0 | 38.13 | 79.63 | 117.76 | ERROR | |
| r-devel-linux-x86_64-debian-gcc | 0.3.0 | 29.18 | 60.56 | 89.74 | ERROR | |
| r-devel-linux-x86_64-fedora-clang | 0.3.0 | 209.14 | OK | |||
| r-devel-linux-x86_64-fedora-gcc | 0.3.0 | 188.77 | OK | |||
| r-devel-windows-ix86+x86_64 | 0.3.0 | 61.00 | 164.00 | 225.00 | OK | |
| r-devel-windows-ix86+x86_64-gcc8 | 0.3.0 | 80.00 | 162.00 | 242.00 | OK | |
| r-patched-linux-x86_64 | 0.3.0 | 30.92 | 110.35 | 141.27 | OK | |
| r-patched-solaris-x86 | 0.3.0 | 270.10 | OK | |||
| r-release-linux-x86_64 | 0.3.0 | 33.79 | 115.20 | 148.99 | OK | |
| r-release-windows-ix86+x86_64 | 0.3.0 | 90.00 | 151.00 | 241.00 | OK | |
| r-release-osx-x86_64 | 0.3.0 | OK | ||||
| r-oldrel-windows-ix86+x86_64 | 0.3.0 | 56.00 | 156.00 | 212.00 | OK | |
| r-oldrel-osx-x86_64 | 0.3.0 | OK |
Version: 0.3.0
Check: examples
Result: ERROR
Running examples in 'BLPestimatoR-Ex.R' failed
The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: estimateBLP
> ### Title: Performs a BLP demand estimation.
> ### Aliases: estimateBLP
>
> ### ** Examples
>
> K<-2 #number of random coefficients
> data <- simulate_BLP_dataset(nmkt = 25, nbrn = 20,
+ Xlin = c("price", "x1", "x2", "x3", "x4", "x5"),
+ Xexo = c("x1", "x2", "x3", "x4", "x5"),
+ Xrandom = paste0("x",1:K),instruments = paste0("iv",1:10),
+ true.parameters = list(Xlin.true.except.price = rep(0.2,5),
+ Xlin.true.price = -0.2,
+ Xrandom.true = rep(2,K),
+ instrument.effects = rep(2,10),
+ instrument.Xexo.effects = rep(1,5)),
+ price.endogeneity = list( mean.xi = -2,
+ mean.eita = 0,
+ cov = cbind( c(1,0.7), c(0.7,1))),
+ printlevel = 0, seed = 234234 )
>
>
> model <- as.formula("shares ~ price + x1 + x2 + x3 + x4 + x5 |
+ x1 + x2 + x3 + x4 + x5 |
+ 0+ x1 + x2 |
+ iv1 + iv2 + iv3 + iv4 + iv5 + iv6 + iv7 + iv8 +iv9 +iv10" )
>
> blp_data <- BLP_data(model = model, market_identifier="cdid",
+ product_id = "prod_id",
+ productData = data,
+ integration_method = "MLHS" ,
+ integration_accuracy = 40,
+ integration_seed = 1)
Mean utility (variable name: `delta`) is initialized with 0 because of missing or invalid par_delta argument.
>
> theta_guesses <- matrix(c(0.5,2), nrow=2)
> rownames(theta_guesses) <- c("x1","x2")
> colnames(theta_guesses) <- "unobs_sd"
>
> blp_est <- estimateBLP(blp_data =blp_data,
+ par_theta2 = theta_guesses,
+ extremumCheck = FALSE ,
+ printLevel = 1 )
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
BLPestimatoR
--- call from context ---
estimateBLP(blp_data = blp_data, par_theta2 = theta_guesses,
extremumCheck = FALSE, printLevel = 1)
--- call from argument ---
if (class(W) == "try-error") stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
--- R stacktrace ---
where 1: estimateBLP(blp_data = blp_data, par_theta2 = theta_guesses,
extremumCheck = FALSE, printLevel = 1)
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (blp_data, par_theta2, solver_method = "BFGS", solver_maxit = 10000,
solver_reltol = 1e-06, standardError = "heteroskedastic",
extremumCheck = FALSE, printLevel = 2, ...)
{
call_arguments <- match.call(expand.dots = TRUE)
nobs <- blp_data$parameters$nobs
K <- blp_data$parameters$K
if (class(blp_data) != "blp_data")
stop("Input has wrong class. Call BLP_data() first.")
Z <- blp_data$data$Z
W <- try(solve((t(Z) %*% Z)))
if (class(W) == "try-error")
stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
xzwz <- t(blp_data$data$X_lin) %*% Z %*% W %*% t(Z)
xzwzx <- xzwz %*% blp_data$data$X_lin
invxzwzx <- try(solve(xzwzx))
if (class(invxzwzx) == "try-error")
stop("Problems with singular matrices. This might be caused by (nearly) linear dependent regressors or weak instruments.")
blp_data$data$W <- W
blp_data$data$xzwz <- xzwz
blp_data$data$invxzwzx <- invxzwzx
start_theta2 <- .prepare_theta2(par_theta2, final_col_names_par = c("unobs_sd",
blp_data$parameters$demographic_names), final_row_names_par = colnames(blp_data$data$X_rand),
K = blp_data$parameters$K, M = blp_data$parameters$total_demogr)
cat("blp_data were prepared with the following arguments:\n")
print(blp_data$call_arguments)
if (printLevel > 0) {
cat("Starting a BLP demand estimation with ", blp_data$parameters$nobs,
" observations in ", blp_data$parameters$nmkt, " markets...\n")
cat("[integration::method", blp_data$integration$integration_method,
" integration::amountDraws", blp_data$integration$amountDraws,
"]\n")
cat("[blp::inner_tol", blp_data$parameters$inner_tol,
" blp::inner_maxit", blp_data$parameters$inner_maxit,
"]\n")
cat("[solver::method", solver_method, " solver::maxit",
solver_maxit, " solver::reltol", solver_reltol, "]\n")
}
blp_results <- new.env(parent = emptyenv())
blp_results$deltaOld <- blp_data$data$delta
blp_results$innerItAll <- c()
blp_results$negShares <- FALSE
blp_results$gradient <- rep(NA_real_, start_theta2$total_par)
start_time <- Sys.time()
res <- optim(par = start_theta2$par_theta2, fn = gmm_obj,
gr = gmm_gr, method = solver_method, control = list(reltol = solver_reltol,
maxit = solver_maxit), indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = printLevel,
...)
solverMessage <- if (res$convergence == 0)
"Successful convergence"
else paste("See error code (optim package)", res$convergence)
outer_it_out <- res$counts[1]
end_time <- Sys.time()
time <- end_time - start_time
cat("------------------------------------------ \n")
cat(paste("Solver message:", solverMessage, "\n"))
if (!(solverMessage == "Successful convergence"))
stop("Cannot compute post estimation results due to failed minimization routine.")
cat("------------------------------------------ \n")
cat("Final GMM evaluation at optimal parameters: \n")
innerItAll_out <- blp_results$innerItAll
blp_results$deltaOld <- rep(0, nobs)
finalTmp <- gmm_obj(par_theta2 = res$par, indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = 3)
delta_out <- blp_results$deltaOld
theta_rc_out <- res$par
theta_lin_out <- blp_results$bet
sij_out <- blp_results$sij
local_min_out <- finalTmp
gradient_out <- blp_results$gradient
jacob_out <- blp_results$jacobian
xi_out <- blp_results$xi
names_rc <- kronecker(start_theta2$final_col_names_par, start_theta2$final_row_names_par,
paste, sep = "*")
relevantRcDem_index <- start_theta2$indices[, "row"] + max(start_theta2$indices[,
"row"]) * (start_theta2$indices[, "col"] - 1)
names(theta_rc_out) <- names_rc[relevantRcDem_index]
X_lin <- blp_data$data$X_lin
Z <- blp_data$data$Z
W <- blp_data$data$W
a <- t(cbind(X_lin, jacob_out)) %*% Z
tmpSE <- try(solve(a %*% W %*% t(a)))
lin_len <- dim(X_lin)[2]
valid_SE <- (standardError %in% c("heteroskedastic", "homoskedastic",
"cluster")) && (length(standardError) == 1)
if (!valid_SE) {
message("Invalid standard error option is provided. Switching to heteroskedastic standard errors...")
standardError <- "heteroskedastic"
}
if (class(tmpSE) == "try-error")
stop("Standard errors cannot be computed due to singular matrices.")
if (standardError == "heteroskedastic") {
cat("Using the heteroskedastic asymptotic variance-covariance matrix... \n")
omega <- xi_out^2
b <- as.matrix(t(Z) %*% Diagonal(length(omega), omega) %*%
Z)
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
if (standardError == "homoskedastic") {
cat("Using the homoskedastic asymptotic variance-covariance matrix... \n")
COV <- c((t(xi_out) %*% xi_out)/nobs) * tmpSE
}
if (standardError == "cluster") {
group_structure <- blp_data$data$group_structure
if (any(is.na(group_structure)) || is.null(group_structure))
stop("Valid group_structure is not availalbe in blp_data. Clustered standard errors require a variable that describes the cluster affiliation.")
group_structure <- data.frame(group_structure = group_structure)
group_matrix <- model.matrix(as.Formula("~0+group_structure"),
group_structure)
tmp <- c(xi_out) * group_matrix
omega <- tmp %*% t(tmp)
b <- t(Z) %*% omega %*% Z
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
seLinear_out <- sqrt(diag(COV))[1:lin_len]
seRc_out <- sqrt(diag(COV))[-(1:lin_len)]
WaldStatistic <- t(matrix(theta_rc_out)) %*% solve(COV[-(1:lin_len),
-(1:lin_len)]) %*% matrix(theta_rc_out)
if (extremumCheck) {
hessian <- invisible(hessian(func = gmm_obj, x = res$par,
indices = start_theta2$indices, blp_results = blp_results,
blp_data = blp_data, printLevel = 0))
hessianEig <- eigen(hessian)$values
isMin_out <- sum(hessianEig > 0) == start_theta2$total_par
isMin_out <- if (isMin_out)
"positive"
else "negative"
cat(paste("Extremum Check:", isMin_out))
}
else {
isMin_out <- NA
}
output <- list(theta_rc = theta_rc_out, theta_lin = theta_lin_out,
se_rc = seRc_out, se_linear = seLinear_out, local_min = local_min_out,
gradient = gradient_out, time = time, outer_it = outer_it_out,
inner_it = innerItAll_out, delta = delta_out, xi = xi_out,
`#demogrCoef` = blp_data$parameters$total_demogr, `#nmkt` = blp_data$parameters$nmkt,
`#nobs` = nobs, `#exoCoef` = length(colnames(blp_data$data$X_exg)),
indices = start_theta2$indices, rand_coef_rownames = start_theta2$final_row_names_par,
rand_coef_colnames = start_theta2$final_col_names_par,
drawsRcMktShape = blp_data$integration$draws_mktShape,
drawsDemMktShape = blp_data$integration$dD, weights = blp_data$integration$weights,
sij = sij_out, WaldStatistic = WaldStatistic, IslocalMin = isMin_out,
outerCrit = solver_reltol, innerCrit = blp_data$parameters$inner_tol,
intMethod = blp_data$integration$method, intdraws = length(blp_data$integration$weights),
standardErrorMethod = standardError, call = call_arguments)
class(output) <- "blp_est"
return(output)
}
<bytecode: 0xaba2cf0>
<environment: namespace:BLPestimatoR>
--- function search by body ---
Function estimateBLP in namespace BLPestimatoR has this body.
----------- END OF FAILURE REPORT --------------
Fatal error: the condition has length > 1
Flavor: r-devel-linux-x86_64-debian-clang
Version: 0.3.0
Check: re-building of vignette outputs
Result: WARN
Error(s) in re-building vignettes:
...
--- re-building 'blp_intro.Rmd' using rmarkdown
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
BLPestimatoR
--- call from context ---
estimateBLP(blp_data = cereal_data, par_theta2 = theta_guesses_cereal,
solver_method = "BFGS", solver_maxit = 1000, solver_reltol = 1e-06,
standardError = "heteroskedastic", extremumCheck = FALSE,
printLevel = 1)
--- call from argument ---
if (class(W) == "try-error") stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
--- R stacktrace ---
where 1: estimateBLP(blp_data = cereal_data, par_theta2 = theta_guesses_cereal,
solver_method = "BFGS", solver_maxit = 1000, solver_reltol = 1e-06,
standardError = "heteroskedastic", extremumCheck = FALSE,
printLevel = 1)
where 2: eval(expr, envir, enclos)
where 3: eval(expr, envir, enclos)
where 4: withVisible(eval(expr, envir, enclos))
where 5: withCallingHandlers(withVisible(eval(expr, envir, enclos)), warning = wHandler,
error = eHandler, message = mHandler)
where 6: handle(ev <- withCallingHandlers(withVisible(eval(expr, envir,
enclos)), warning = wHandler, error = eHandler, message = mHandler))
where 7: timing_fn(handle(ev <- withCallingHandlers(withVisible(eval(expr,
envir, enclos)), warning = wHandler, error = eHandler, message = mHandler)))
where 8: evaluate_call(expr, parsed$src[[i]], envir = envir, enclos = enclos,
debug = debug, last = i == length(out), use_try = stop_on_error !=
2L, keep_warning = keep_warning, keep_message = keep_message,
output_handler = output_handler, include_timing = include_timing)
where 9: evaluate::evaluate(...)
where 10: evaluate(code, envir = env, new_device = FALSE, keep_warning = !isFALSE(options$warning),
keep_message = !isFALSE(options$message), stop_on_error = if (options$error &&
options$include) 0L else 2L, output_handler = knit_handlers(options$render,
options))
where 11: in_dir(input_dir(), evaluate(code, envir = env, new_device = FALSE,
keep_warning = !isFALSE(options$warning), keep_message = !isFALSE(options$message),
stop_on_error = if (options$error && options$include) 0L else 2L,
output_handler = knit_handlers(options$render, options)))
where 12: block_exec(params)
where 13: call_block(x)
where 14: process_group.block(group)
where 15: process_group(group)
where 16: withCallingHandlers(if (tangle) process_tangle(group) else process_group(group),
error = function(e) {
setwd(wd)
cat(res, sep = "\n", file = output %n% "")
message("Quitting from lines ", paste(current_lines(i),
collapse = "-"), " (", knit_concord$get("infile"),
") ")
})
where 17: process_file(text, output)
where 18: knitr::knit(knit_input, knit_output, envir = envir, quiet = quiet,
encoding = encoding)
where 19: rmarkdown::render(file, encoding = encoding, quiet = quiet, envir = globalenv(),
...)
where 20: vweave_rmarkdown(...)
where 21: engine$weave(file, quiet = quiet, encoding = enc)
where 22: doTryCatch(return(expr), name, parentenv, handler)
where 23: tryCatchOne(expr, names, parentenv, handlers[[1L]])
where 24: tryCatchList(expr, classes, parentenv, handlers)
where 25: tryCatch({
engine$weave(file, quiet = quiet, encoding = enc)
setwd(startdir)
output <- find_vignette_product(name, by = "weave", engine = engine)
if (!have.makefile && vignette_is_tex(output)) {
texi2pdf(file = output, clean = FALSE, quiet = quiet)
output <- find_vignette_product(name, by = "texi2pdf",
engine = engine)
}
outputs <- c(outputs, output)
}, error = function(e) {
thisOK <<- FALSE
fails <<- c(fails, file)
message(gettextf("Error: processing vignette '%s' failed with diagnostics:\n%s",
file, conditionMessage(e)))
})
where 26: tools:::buildVignettes(dir = "/home/hornik/tmp/R.check/r-devel-clang/Work/PKGS/BLPestimatoR.Rcheck/vign_test/BLPestimatoR",
ser_elibs = "/tmp/RtmpS7rUfb/file36e5885ac7e.rds")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (blp_data, par_theta2, solver_method = "BFGS", solver_maxit = 10000,
solver_reltol = 1e-06, standardError = "heteroskedastic",
extremumCheck = FALSE, printLevel = 2, ...)
{
call_arguments <- match.call(expand.dots = TRUE)
nobs <- blp_data$parameters$nobs
K <- blp_data$parameters$K
if (class(blp_data) != "blp_data")
stop("Input has wrong class. Call BLP_data() first.")
Z <- blp_data$data$Z
W <- try(solve((t(Z) %*% Z)))
if (class(W) == "try-error")
stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
xzwz <- t(blp_data$data$X_lin) %*% Z %*% W %*% t(Z)
xzwzx <- xzwz %*% blp_data$data$X_lin
invxzwzx <- try(solve(xzwzx))
if (class(invxzwzx) == "try-error")
stop("Problems with singular matrices. This might be caused by (nearly) linear dependent regressors or weak instruments.")
blp_data$data$W <- W
blp_data$data$xzwz <- xzwz
blp_data$data$invxzwzx <- invxzwzx
start_theta2 <- .prepare_theta2(par_theta2, final_col_names_par = c("unobs_sd",
blp_data$parameters$demographic_names), final_row_names_par = colnames(blp_data$data$X_rand),
K = blp_data$parameters$K, M = blp_data$parameters$total_demogr)
cat("blp_data were prepared with the following arguments:\n")
print(blp_data$call_arguments)
if (printLevel > 0) {
cat("Starting a BLP demand estimation with ", blp_data$parameters$nobs,
" observations in ", blp_data$parameters$nmkt, " markets...\n")
cat("[integration::method", blp_data$integration$integration_method,
" integration::amountDraws", blp_data$integration$amountDraws,
"]\n")
cat("[blp::inner_tol", blp_data$parameters$inner_tol,
" blp::inner_maxit", blp_data$parameters$inner_maxit,
"]\n")
cat("[solver::method", solver_method, " solver::maxit",
solver_maxit, " solver::reltol", solver_reltol, "]\n")
}
blp_results <- new.env(parent = emptyenv())
blp_results$deltaOld <- blp_data$data$delta
blp_results$innerItAll <- c()
blp_results$negShares <- FALSE
blp_results$gradient <- rep(NA_real_, start_theta2$total_par)
start_time <- Sys.time()
res <- optim(par = start_theta2$par_theta2, fn = gmm_obj,
gr = gmm_gr, method = solver_method, control = list(reltol = solver_reltol,
maxit = solver_maxit), indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = printLevel,
...)
solverMessage <- if (res$convergence == 0)
"Successful convergence"
else paste("See error code (optim package)", res$convergence)
outer_it_out <- res$counts[1]
end_time <- Sys.time()
time <- end_time - start_time
cat("------------------------------------------ \n")
cat(paste("Solver message:", solverMessage, "\n"))
if (!(solverMessage == "Successful convergence"))
stop("Cannot compute post estimation results due to failed minimization routine.")
cat("------------------------------------------ \n")
cat("Final GMM evaluation at optimal parameters: \n")
innerItAll_out <- blp_results$innerItAll
blp_results$deltaOld <- rep(0, nobs)
finalTmp <- gmm_obj(par_theta2 = res$par, indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = 3)
delta_out <- blp_results$deltaOld
theta_rc_out <- res$par
theta_lin_out <- blp_results$bet
sij_out <- blp_results$sij
local_min_out <- finalTmp
gradient_out <- blp_results$gradient
jacob_out <- blp_results$jacobian
xi_out <- blp_results$xi
names_rc <- kronecker(start_theta2$final_col_names_par, start_theta2$final_row_names_par,
paste, sep = "*")
relevantRcDem_index <- start_theta2$indices[, "row"] + max(start_theta2$indices[,
"row"]) * (start_theta2$indices[, "col"] - 1)
names(theta_rc_out) <- names_rc[relevantRcDem_index]
X_lin <- blp_data$data$X_lin
Z <- blp_data$data$Z
W <- blp_data$data$W
a <- t(cbind(X_lin, jacob_out)) %*% Z
tmpSE <- try(solve(a %*% W %*% t(a)))
lin_len <- dim(X_lin)[2]
valid_SE <- (standardError %in% c("heteroskedastic", "homoskedastic",
"cluster")) && (length(standardError) == 1)
if (!valid_SE) {
message("Invalid standard error option is provided. Switching to heteroskedastic standard errors...")
standardError <- "heteroskedastic"
}
if (class(tmpSE) == "try-error")
stop("Standard errors cannot be computed due to singular matrices.")
if (standardError == "heteroskedastic") {
cat("Using the heteroskedastic asymptotic variance-covariance matrix... \n")
omega <- xi_out^2
b <- as.matrix(t(Z) %*% Diagonal(length(omega), omega) %*%
Z)
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
if (standardError == "homoskedastic") {
cat("Using the homoskedastic asymptotic variance-covariance matrix... \n")
COV <- c((t(xi_out) %*% xi_out)/nobs) * tmpSE
}
if (standardError == "cluster") {
group_structure <- blp_data$data$group_structure
if (any(is.na(group_structure)) || is.null(group_structure))
stop("Valid group_structure is not availalbe in blp_data. Clustered standard errors require a variable that describes the cluster affiliation.")
group_structure <- data.frame(group_structure = group_structure)
group_matrix <- model.matrix(as.Formula("~0+group_structure"),
group_structure)
tmp <- c(xi_out) * group_matrix
omega <- tmp %*% t(tmp)
b <- t(Z) %*% omega %*% Z
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
seLinear_out <- sqrt(diag(COV))[1:lin_len]
seRc_out <- sqrt(diag(COV))[-(1:lin_len)]
WaldStatistic <- t(matrix(theta_rc_out)) %*% solve(COV[-(1:lin_len),
-(1:lin_len)]) %*% matrix(theta_rc_out)
if (extremumCheck) {
hessian <- invisible(hessian(func = gmm_obj, x = res$par,
indices = start_theta2$indices, blp_results = blp_results,
blp_data = blp_data, printLevel = 0))
hessianEig <- eigen(hessian)$values
isMin_out <- sum(hessianEig > 0) == start_theta2$total_par
isMin_out <- if (isMin_out)
"positive"
else "negative"
cat(paste("Extremum Check:", isMin_out))
}
else {
isMin_out <- NA
}
output <- list(theta_rc = theta_rc_out, theta_lin = theta_lin_out,
se_rc = seRc_out, se_linear = seLinear_out, local_min = local_min_out,
gradient = gradient_out, time = time, outer_it = outer_it_out,
inner_it = innerItAll_out, delta = delta_out, xi = xi_out,
`#demogrCoef` = blp_data$parameters$total_demogr, `#nmkt` = blp_data$parameters$nmkt,
`#nobs` = nobs, `#exoCoef` = length(colnames(blp_data$data$X_exg)),
indices = start_theta2$indices, rand_coef_rownames = start_theta2$final_row_names_par,
rand_coef_colnames = start_theta2$final_col_names_par,
drawsRcMktShape = blp_data$integration$draws_mktShape,
drawsDemMktShape = blp_data$integration$dD, weights = blp_data$integration$weights,
sij = sij_out, WaldStatistic = WaldStatistic, IslocalMin = isMin_out,
outerCrit = solver_reltol, innerCrit = blp_data$parameters$inner_tol,
intMethod = blp_data$integration$method, intdraws = length(blp_data$integration$weights),
standardErrorMethod = standardError, call = call_arguments)
class(output) <- "blp_est"
return(output)
}
<bytecode: 0x5573c18>
<environment: namespace:BLPestimatoR>
--- function search by body ---
Function estimateBLP in namespace BLPestimatoR has this body.
----------- END OF FAILURE REPORT --------------
Fatal error: the condition has length > 1
Flavor: r-devel-linux-x86_64-debian-clang
Version: 0.3.0
Check: examples
Result: ERROR
Running examples in ‘BLPestimatoR-Ex.R’ failed
The error most likely occurred in:
> base::assign(".ptime", proc.time(), pos = "CheckExEnv")
> ### Name: estimateBLP
> ### Title: Performs a BLP demand estimation.
> ### Aliases: estimateBLP
>
> ### ** Examples
>
> K<-2 #number of random coefficients
> data <- simulate_BLP_dataset(nmkt = 25, nbrn = 20,
+ Xlin = c("price", "x1", "x2", "x3", "x4", "x5"),
+ Xexo = c("x1", "x2", "x3", "x4", "x5"),
+ Xrandom = paste0("x",1:K),instruments = paste0("iv",1:10),
+ true.parameters = list(Xlin.true.except.price = rep(0.2,5),
+ Xlin.true.price = -0.2,
+ Xrandom.true = rep(2,K),
+ instrument.effects = rep(2,10),
+ instrument.Xexo.effects = rep(1,5)),
+ price.endogeneity = list( mean.xi = -2,
+ mean.eita = 0,
+ cov = cbind( c(1,0.7), c(0.7,1))),
+ printlevel = 0, seed = 234234 )
>
>
> model <- as.formula("shares ~ price + x1 + x2 + x3 + x4 + x5 |
+ x1 + x2 + x3 + x4 + x5 |
+ 0+ x1 + x2 |
+ iv1 + iv2 + iv3 + iv4 + iv5 + iv6 + iv7 + iv8 +iv9 +iv10" )
>
> blp_data <- BLP_data(model = model, market_identifier="cdid",
+ product_id = "prod_id",
+ productData = data,
+ integration_method = "MLHS" ,
+ integration_accuracy = 40,
+ integration_seed = 1)
Mean utility (variable name: `delta`) is initialized with 0 because of missing or invalid par_delta argument.
>
> theta_guesses <- matrix(c(0.5,2), nrow=2)
> rownames(theta_guesses) <- c("x1","x2")
> colnames(theta_guesses) <- "unobs_sd"
>
> blp_est <- estimateBLP(blp_data =blp_data,
+ par_theta2 = theta_guesses,
+ extremumCheck = FALSE ,
+ printLevel = 1 )
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
BLPestimatoR
--- call from context ---
estimateBLP(blp_data = blp_data, par_theta2 = theta_guesses,
extremumCheck = FALSE, printLevel = 1)
--- call from argument ---
if (class(W) == "try-error") stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
--- R stacktrace ---
where 1: estimateBLP(blp_data = blp_data, par_theta2 = theta_guesses,
extremumCheck = FALSE, printLevel = 1)
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (blp_data, par_theta2, solver_method = "BFGS", solver_maxit = 10000,
solver_reltol = 1e-06, standardError = "heteroskedastic",
extremumCheck = FALSE, printLevel = 2, ...)
{
call_arguments <- match.call(expand.dots = TRUE)
nobs <- blp_data$parameters$nobs
K <- blp_data$parameters$K
if (class(blp_data) != "blp_data")
stop("Input has wrong class. Call BLP_data() first.")
Z <- blp_data$data$Z
W <- try(solve((t(Z) %*% Z)))
if (class(W) == "try-error")
stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
xzwz <- t(blp_data$data$X_lin) %*% Z %*% W %*% t(Z)
xzwzx <- xzwz %*% blp_data$data$X_lin
invxzwzx <- try(solve(xzwzx))
if (class(invxzwzx) == "try-error")
stop("Problems with singular matrices. This might be caused by (nearly) linear dependent regressors or weak instruments.")
blp_data$data$W <- W
blp_data$data$xzwz <- xzwz
blp_data$data$invxzwzx <- invxzwzx
start_theta2 <- .prepare_theta2(par_theta2, final_col_names_par = c("unobs_sd",
blp_data$parameters$demographic_names), final_row_names_par = colnames(blp_data$data$X_rand),
K = blp_data$parameters$K, M = blp_data$parameters$total_demogr)
cat("blp_data were prepared with the following arguments:\n")
print(blp_data$call_arguments)
if (printLevel > 0) {
cat("Starting a BLP demand estimation with ", blp_data$parameters$nobs,
" observations in ", blp_data$parameters$nmkt, " markets...\n")
cat("[integration::method", blp_data$integration$integration_method,
" integration::amountDraws", blp_data$integration$amountDraws,
"]\n")
cat("[blp::inner_tol", blp_data$parameters$inner_tol,
" blp::inner_maxit", blp_data$parameters$inner_maxit,
"]\n")
cat("[solver::method", solver_method, " solver::maxit",
solver_maxit, " solver::reltol", solver_reltol, "]\n")
}
blp_results <- new.env(parent = emptyenv())
blp_results$deltaOld <- blp_data$data$delta
blp_results$innerItAll <- c()
blp_results$negShares <- FALSE
blp_results$gradient <- rep(NA_real_, start_theta2$total_par)
start_time <- Sys.time()
res <- optim(par = start_theta2$par_theta2, fn = gmm_obj,
gr = gmm_gr, method = solver_method, control = list(reltol = solver_reltol,
maxit = solver_maxit), indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = printLevel,
...)
solverMessage <- if (res$convergence == 0)
"Successful convergence"
else paste("See error code (optim package)", res$convergence)
outer_it_out <- res$counts[1]
end_time <- Sys.time()
time <- end_time - start_time
cat("------------------------------------------ \n")
cat(paste("Solver message:", solverMessage, "\n"))
if (!(solverMessage == "Successful convergence"))
stop("Cannot compute post estimation results due to failed minimization routine.")
cat("------------------------------------------ \n")
cat("Final GMM evaluation at optimal parameters: \n")
innerItAll_out <- blp_results$innerItAll
blp_results$deltaOld <- rep(0, nobs)
finalTmp <- gmm_obj(par_theta2 = res$par, indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = 3)
delta_out <- blp_results$deltaOld
theta_rc_out <- res$par
theta_lin_out <- blp_results$bet
sij_out <- blp_results$sij
local_min_out <- finalTmp
gradient_out <- blp_results$gradient
jacob_out <- blp_results$jacobian
xi_out <- blp_results$xi
names_rc <- kronecker(start_theta2$final_col_names_par, start_theta2$final_row_names_par,
paste, sep = "*")
relevantRcDem_index <- start_theta2$indices[, "row"] + max(start_theta2$indices[,
"row"]) * (start_theta2$indices[, "col"] - 1)
names(theta_rc_out) <- names_rc[relevantRcDem_index]
X_lin <- blp_data$data$X_lin
Z <- blp_data$data$Z
W <- blp_data$data$W
a <- t(cbind(X_lin, jacob_out)) %*% Z
tmpSE <- try(solve(a %*% W %*% t(a)))
lin_len <- dim(X_lin)[2]
valid_SE <- (standardError %in% c("heteroskedastic", "homoskedastic",
"cluster")) && (length(standardError) == 1)
if (!valid_SE) {
message("Invalid standard error option is provided. Switching to heteroskedastic standard errors...")
standardError <- "heteroskedastic"
}
if (class(tmpSE) == "try-error")
stop("Standard errors cannot be computed due to singular matrices.")
if (standardError == "heteroskedastic") {
cat("Using the heteroskedastic asymptotic variance-covariance matrix... \n")
omega <- xi_out^2
b <- as.matrix(t(Z) %*% Diagonal(length(omega), omega) %*%
Z)
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
if (standardError == "homoskedastic") {
cat("Using the homoskedastic asymptotic variance-covariance matrix... \n")
COV <- c((t(xi_out) %*% xi_out)/nobs) * tmpSE
}
if (standardError == "cluster") {
group_structure <- blp_data$data$group_structure
if (any(is.na(group_structure)) || is.null(group_structure))
stop("Valid group_structure is not availalbe in blp_data. Clustered standard errors require a variable that describes the cluster affiliation.")
group_structure <- data.frame(group_structure = group_structure)
group_matrix <- model.matrix(as.Formula("~0+group_structure"),
group_structure)
tmp <- c(xi_out) * group_matrix
omega <- tmp %*% t(tmp)
b <- t(Z) %*% omega %*% Z
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
seLinear_out <- sqrt(diag(COV))[1:lin_len]
seRc_out <- sqrt(diag(COV))[-(1:lin_len)]
WaldStatistic <- t(matrix(theta_rc_out)) %*% solve(COV[-(1:lin_len),
-(1:lin_len)]) %*% matrix(theta_rc_out)
if (extremumCheck) {
hessian <- invisible(hessian(func = gmm_obj, x = res$par,
indices = start_theta2$indices, blp_results = blp_results,
blp_data = blp_data, printLevel = 0))
hessianEig <- eigen(hessian)$values
isMin_out <- sum(hessianEig > 0) == start_theta2$total_par
isMin_out <- if (isMin_out)
"positive"
else "negative"
cat(paste("Extremum Check:", isMin_out))
}
else {
isMin_out <- NA
}
output <- list(theta_rc = theta_rc_out, theta_lin = theta_lin_out,
se_rc = seRc_out, se_linear = seLinear_out, local_min = local_min_out,
gradient = gradient_out, time = time, outer_it = outer_it_out,
inner_it = innerItAll_out, delta = delta_out, xi = xi_out,
`#demogrCoef` = blp_data$parameters$total_demogr, `#nmkt` = blp_data$parameters$nmkt,
`#nobs` = nobs, `#exoCoef` = length(colnames(blp_data$data$X_exg)),
indices = start_theta2$indices, rand_coef_rownames = start_theta2$final_row_names_par,
rand_coef_colnames = start_theta2$final_col_names_par,
drawsRcMktShape = blp_data$integration$draws_mktShape,
drawsDemMktShape = blp_data$integration$dD, weights = blp_data$integration$weights,
sij = sij_out, WaldStatistic = WaldStatistic, IslocalMin = isMin_out,
outerCrit = solver_reltol, innerCrit = blp_data$parameters$inner_tol,
intMethod = blp_data$integration$method, intdraws = length(blp_data$integration$weights),
standardErrorMethod = standardError, call = call_arguments)
class(output) <- "blp_est"
return(output)
}
<bytecode: 0x55aeb965d168>
<environment: namespace:BLPestimatoR>
--- function search by body ---
Function estimateBLP in namespace BLPestimatoR has this body.
----------- END OF FAILURE REPORT --------------
Fatal error: the condition has length > 1
Flavor: r-devel-linux-x86_64-debian-gcc
Version: 0.3.0
Check: re-building of vignette outputs
Result: WARN
Error(s) in re-building vignettes:
...
--- re-building ‘blp_intro.Rmd’ using rmarkdown
----------- FAILURE REPORT --------------
--- failure: the condition has length > 1 ---
--- srcref ---
:
--- package (from environment) ---
BLPestimatoR
--- call from context ---
estimateBLP(blp_data = cereal_data, par_theta2 = theta_guesses_cereal,
solver_method = "BFGS", solver_maxit = 1000, solver_reltol = 1e-06,
standardError = "heteroskedastic", extremumCheck = FALSE,
printLevel = 1)
--- call from argument ---
if (class(W) == "try-error") stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
--- R stacktrace ---
where 1: estimateBLP(blp_data = cereal_data, par_theta2 = theta_guesses_cereal,
solver_method = "BFGS", solver_maxit = 1000, solver_reltol = 1e-06,
standardError = "heteroskedastic", extremumCheck = FALSE,
printLevel = 1)
where 2: eval(expr, envir, enclos)
where 3: eval(expr, envir, enclos)
where 4: withVisible(eval(expr, envir, enclos))
where 5: withCallingHandlers(withVisible(eval(expr, envir, enclos)), warning = wHandler,
error = eHandler, message = mHandler)
where 6: handle(ev <- withCallingHandlers(withVisible(eval(expr, envir,
enclos)), warning = wHandler, error = eHandler, message = mHandler))
where 7: timing_fn(handle(ev <- withCallingHandlers(withVisible(eval(expr,
envir, enclos)), warning = wHandler, error = eHandler, message = mHandler)))
where 8: evaluate_call(expr, parsed$src[[i]], envir = envir, enclos = enclos,
debug = debug, last = i == length(out), use_try = stop_on_error !=
2L, keep_warning = keep_warning, keep_message = keep_message,
output_handler = output_handler, include_timing = include_timing)
where 9: evaluate::evaluate(...)
where 10: evaluate(code, envir = env, new_device = FALSE, keep_warning = !isFALSE(options$warning),
keep_message = !isFALSE(options$message), stop_on_error = if (options$error &&
options$include) 0L else 2L, output_handler = knit_handlers(options$render,
options))
where 11: in_dir(input_dir(), evaluate(code, envir = env, new_device = FALSE,
keep_warning = !isFALSE(options$warning), keep_message = !isFALSE(options$message),
stop_on_error = if (options$error && options$include) 0L else 2L,
output_handler = knit_handlers(options$render, options)))
where 12: block_exec(params)
where 13: call_block(x)
where 14: process_group.block(group)
where 15: process_group(group)
where 16: withCallingHandlers(if (tangle) process_tangle(group) else process_group(group),
error = function(e) {
setwd(wd)
cat(res, sep = "\n", file = output %n% "")
message("Quitting from lines ", paste(current_lines(i),
collapse = "-"), " (", knit_concord$get("infile"),
") ")
})
where 17: process_file(text, output)
where 18: knitr::knit(knit_input, knit_output, envir = envir, quiet = quiet,
encoding = encoding)
where 19: rmarkdown::render(file, encoding = encoding, quiet = quiet, envir = globalenv(),
...)
where 20: vweave_rmarkdown(...)
where 21: engine$weave(file, quiet = quiet, encoding = enc)
where 22: doTryCatch(return(expr), name, parentenv, handler)
where 23: tryCatchOne(expr, names, parentenv, handlers[[1L]])
where 24: tryCatchList(expr, classes, parentenv, handlers)
where 25: tryCatch({
engine$weave(file, quiet = quiet, encoding = enc)
setwd(startdir)
output <- find_vignette_product(name, by = "weave", engine = engine)
if (!have.makefile && vignette_is_tex(output)) {
texi2pdf(file = output, clean = FALSE, quiet = quiet)
output <- find_vignette_product(name, by = "texi2pdf",
engine = engine)
}
outputs <- c(outputs, output)
}, error = function(e) {
thisOK <<- FALSE
fails <<- c(fails, file)
message(gettextf("Error: processing vignette '%s' failed with diagnostics:\n%s",
file, conditionMessage(e)))
})
where 26: tools:::buildVignettes(dir = "/home/hornik/tmp/R.check/r-devel-gcc/Work/PKGS/BLPestimatoR.Rcheck/vign_test/BLPestimatoR",
ser_elibs = "/home/hornik/tmp/scratch/Rtmpct4jBw/file2c6b56bcbdd3.rds")
--- value of length: 2 type: logical ---
[1] FALSE FALSE
--- function from context ---
function (blp_data, par_theta2, solver_method = "BFGS", solver_maxit = 10000,
solver_reltol = 1e-06, standardError = "heteroskedastic",
extremumCheck = FALSE, printLevel = 2, ...)
{
call_arguments <- match.call(expand.dots = TRUE)
nobs <- blp_data$parameters$nobs
K <- blp_data$parameters$K
if (class(blp_data) != "blp_data")
stop("Input has wrong class. Call BLP_data() first.")
Z <- blp_data$data$Z
W <- try(solve((t(Z) %*% Z)))
if (class(W) == "try-error")
stop("Problems with singular matrizes. This might be caused by (nearly) linear dependent regressors or weak instruments.")
xzwz <- t(blp_data$data$X_lin) %*% Z %*% W %*% t(Z)
xzwzx <- xzwz %*% blp_data$data$X_lin
invxzwzx <- try(solve(xzwzx))
if (class(invxzwzx) == "try-error")
stop("Problems with singular matrices. This might be caused by (nearly) linear dependent regressors or weak instruments.")
blp_data$data$W <- W
blp_data$data$xzwz <- xzwz
blp_data$data$invxzwzx <- invxzwzx
start_theta2 <- .prepare_theta2(par_theta2, final_col_names_par = c("unobs_sd",
blp_data$parameters$demographic_names), final_row_names_par = colnames(blp_data$data$X_rand),
K = blp_data$parameters$K, M = blp_data$parameters$total_demogr)
cat("blp_data were prepared with the following arguments:\n")
print(blp_data$call_arguments)
if (printLevel > 0) {
cat("Starting a BLP demand estimation with ", blp_data$parameters$nobs,
" observations in ", blp_data$parameters$nmkt, " markets...\n")
cat("[integration::method", blp_data$integration$integration_method,
" integration::amountDraws", blp_data$integration$amountDraws,
"]\n")
cat("[blp::inner_tol", blp_data$parameters$inner_tol,
" blp::inner_maxit", blp_data$parameters$inner_maxit,
"]\n")
cat("[solver::method", solver_method, " solver::maxit",
solver_maxit, " solver::reltol", solver_reltol, "]\n")
}
blp_results <- new.env(parent = emptyenv())
blp_results$deltaOld <- blp_data$data$delta
blp_results$innerItAll <- c()
blp_results$negShares <- FALSE
blp_results$gradient <- rep(NA_real_, start_theta2$total_par)
start_time <- Sys.time()
res <- optim(par = start_theta2$par_theta2, fn = gmm_obj,
gr = gmm_gr, method = solver_method, control = list(reltol = solver_reltol,
maxit = solver_maxit), indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = printLevel,
...)
solverMessage <- if (res$convergence == 0)
"Successful convergence"
else paste("See error code (optim package)", res$convergence)
outer_it_out <- res$counts[1]
end_time <- Sys.time()
time <- end_time - start_time
cat("------------------------------------------ \n")
cat(paste("Solver message:", solverMessage, "\n"))
if (!(solverMessage == "Successful convergence"))
stop("Cannot compute post estimation results due to failed minimization routine.")
cat("------------------------------------------ \n")
cat("Final GMM evaluation at optimal parameters: \n")
innerItAll_out <- blp_results$innerItAll
blp_results$deltaOld <- rep(0, nobs)
finalTmp <- gmm_obj(par_theta2 = res$par, indices = start_theta2$indices,
blp_results = blp_results, blp_data = blp_data, printLevel = 3)
delta_out <- blp_results$deltaOld
theta_rc_out <- res$par
theta_lin_out <- blp_results$bet
sij_out <- blp_results$sij
local_min_out <- finalTmp
gradient_out <- blp_results$gradient
jacob_out <- blp_results$jacobian
xi_out <- blp_results$xi
names_rc <- kronecker(start_theta2$final_col_names_par, start_theta2$final_row_names_par,
paste, sep = "*")
relevantRcDem_index <- start_theta2$indices[, "row"] + max(start_theta2$indices[,
"row"]) * (start_theta2$indices[, "col"] - 1)
names(theta_rc_out) <- names_rc[relevantRcDem_index]
X_lin <- blp_data$data$X_lin
Z <- blp_data$data$Z
W <- blp_data$data$W
a <- t(cbind(X_lin, jacob_out)) %*% Z
tmpSE <- try(solve(a %*% W %*% t(a)))
lin_len <- dim(X_lin)[2]
valid_SE <- (standardError %in% c("heteroskedastic", "homoskedastic",
"cluster")) && (length(standardError) == 1)
if (!valid_SE) {
message("Invalid standard error option is provided. Switching to heteroskedastic standard errors...")
standardError <- "heteroskedastic"
}
if (class(tmpSE) == "try-error")
stop("Standard errors cannot be computed due to singular matrices.")
if (standardError == "heteroskedastic") {
cat("Using the heteroskedastic asymptotic variance-covariance matrix... \n")
omega <- xi_out^2
b <- as.matrix(t(Z) %*% Diagonal(length(omega), omega) %*%
Z)
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
if (standardError == "homoskedastic") {
cat("Using the homoskedastic asymptotic variance-covariance matrix... \n")
COV <- c((t(xi_out) %*% xi_out)/nobs) * tmpSE
}
if (standardError == "cluster") {
group_structure <- blp_data$data$group_structure
if (any(is.na(group_structure)) || is.null(group_structure))
stop("Valid group_structure is not availalbe in blp_data. Clustered standard errors require a variable that describes the cluster affiliation.")
group_structure <- data.frame(group_structure = group_structure)
group_matrix <- model.matrix(as.Formula("~0+group_structure"),
group_structure)
tmp <- c(xi_out) * group_matrix
omega <- tmp %*% t(tmp)
b <- t(Z) %*% omega %*% Z
COV <- tmpSE %*% a %*% W %*% b %*% W %*% t(a) %*% tmpSE
}
seLinear_out <- sqrt(diag(COV))[1:lin_len]
seRc_out <- sqrt(diag(COV))[-(1:lin_len)]
WaldStatistic <- t(matrix(theta_rc_out)) %*% solve(COV[-(1:lin_len),
-(1:lin_len)]) %*% matrix(theta_rc_out)
if (extremumCheck) {
hessian <- invisible(hessian(func = gmm_obj, x = res$par,
indices = start_theta2$indices, blp_results = blp_results,
blp_data = blp_data, printLevel = 0))
hessianEig <- eigen(hessian)$values
isMin_out <- sum(hessianEig > 0) == start_theta2$total_par
isMin_out <- if (isMin_out)
"positive"
else "negative"
cat(paste("Extremum Check:", isMin_out))
}
else {
isMin_out <- NA
}
output <- list(theta_rc = theta_rc_out, theta_lin = theta_lin_out,
se_rc = seRc_out, se_linear = seLinear_out, local_min = local_min_out,
gradient = gradient_out, time = time, outer_it = outer_it_out,
inner_it = innerItAll_out, delta = delta_out, xi = xi_out,
`#demogrCoef` = blp_data$parameters$total_demogr, `#nmkt` = blp_data$parameters$nmkt,
`#nobs` = nobs, `#exoCoef` = length(colnames(blp_data$data$X_exg)),
indices = start_theta2$indices, rand_coef_rownames = start_theta2$final_row_names_par,
rand_coef_colnames = start_theta2$final_col_names_par,
drawsRcMktShape = blp_data$integration$draws_mktShape,
drawsDemMktShape = blp_data$integration$dD, weights = blp_data$integration$weights,
sij = sij_out, WaldStatistic = WaldStatistic, IslocalMin = isMin_out,
outerCrit = solver_reltol, innerCrit = blp_data$parameters$inner_tol,
intMethod = blp_data$integration$method, intdraws = length(blp_data$integration$weights),
standardErrorMethod = standardError, call = call_arguments)
class(output) <- "blp_est"
return(output)
}
<bytecode: 0x561ea8601160>
<environment: namespace:BLPestimatoR>
--- function search by body ---
Function estimateBLP in namespace BLPestimatoR has this body.
----------- END OF FAILURE REPORT --------------
Fatal error: the condition has length > 1
Flavor: r-devel-linux-x86_64-debian-gcc