## ----include = FALSE----------------------------------------------------------
knitr::opts_chunk$set(
  collapse = TRUE,
  comment = "#>"
)

## ----setup, include=FALSE-----------------------------------------------------
library(fastqrs)

## ----load-libraries, echo=TRUE------------------------------------------------
# Load the package
library("sampleSelection")
library("quantreg")
library("fastqrs")
library("copula")
library("ggplot2")

## ----load-data, echo=TRUE-----------------------------------------------------
data("Mroz87")
Mroz87$kids <- ( Mroz87$kids5 + Mroz87$kids618 > 0 )

y <- Mroz87$wage
d <- Mroz87$lfp
x <- cbind(Mroz87$exper, Mroz87$exper^2,Mroz87$educ,
           Mroz87$city)
z <- cbind(Mroz87$age, Mroz87$age^2, Mroz87$faminc,
           Mroz87$kids, Mroz87$educ)
w <- rep(1,NROW(y))

## ----set-parameters, echo=TRUE------------------------------------------------
link <- "probit"
family <- "Gaussian"
gridtheta <- seq(from = -1, to = 1, by = .05)
Q1 <- 9
Q2 <- 19
P <- 10
m <- 1

## ----estimation, echo=TRUE----------------------------------------------------
qrs <- qrs.fast(y, x, d, z, w, Q1, Q2, P, link, family,
                gridtheta, m)

## ----graph1, echo=TRUE, fig=TRUE, fig.width=6, fig.height=4-------------------
data_for_plot <- data.frame(
   Quantile = seq(1/(Q2+1), Q2/(Q2+1), length.out = Q2),
   Coefficient = as.numeric(qrs$beta[5, ]))
ggplot(data_for_plot, aes(x = Quantile, y = Coefficient)) +
   geom_line(color = "blue", linewidth = 1) +
   geom_point(color = "red", size = 2) +
   labs(
       title = "QRS Coefficients for Urban",
       x = "Quantile",
       y = "Estimated Coefficient"
   )

## ----bootstrap, echo=TRUE-----------------------------------------------------
if (interactive()) {
  set.seed(1)
  reps<-200
  alpha <- c(0.1,0.05,0.01)
  qrs.bt <- qrs.fast.bt(y, x, d, z, w, Q1, Q2, P, link, family,
                        gridtheta, m, qrs$b1, reps, alpha)
} else {
  #bootstrap_results <- load(system.file("extdata", "bootstrap_results.RData", package="fastqrs"))
  bootstrap_file <- system.file("extdata", "bootstrap_results.RData", package="fastqrs")
  if (bootstrap_file != "") {
    load(bootstrap_file)
  } else {
    stop("Bootstrap results file not found. Make sure the package is correctly installed.")
  }
}

## ----graph2, echo=TRUE, fig=TRUE, fig.width=6, fig.height=4-------------------
data_for_plot <- data.frame(
   Quantile = seq(1/(Q2+1), Q2/(Q2+1), length.out = Q2),
   Coefficient = as.numeric(qrs$beta[5, ]),
   Upper = as.numeric(qrs.bt$betaub[5, ,1]),
   Lower = as.numeric(qrs.bt$betalb[5, ,1])
)
ggplot(data_for_plot, aes(x = Quantile)) +
   # Line for coefficients
   geom_line(aes(y = Coefficient), color = "blue", linewidth = 1) +
   geom_point(aes(y = Coefficient), color = "blue", size = 2) +
   # Lines for confidence bands
   geom_line(aes(y = Upper), color = "red",
             linetype = "dashed", linewidth = 0.8) +
   geom_line(aes(y = Lower), color = "red",
             linetype = "dashed", linewidth = 0.8) +
   # Add labels and title
   labs(
   title = "QRS Coefficients for Urban with 90% Confidence Bands",
   x = "Quantile",
   y = "Estimated Coefficient"
   )