install.packages("copula")
install.packages("moments")

library(copula)
library(moments)

set.seed(123)

n <- 1000
m <- 2000

qWeib_custom <- function(u, xi, k) (-log1p(-u) / xi)^(1 / k)

rBWRR2 <- function(n, alpha, k, xi1, xi2, b1, b2) {
  stopifnot(alpha > 0, alpha <= 1, k > 0, xi1 > 0, xi2 > 0, b1 > 0, b2 > 0)
  cop <- gumbelCopula(param = 1/alpha, dim = 2)
  UV  <- rCopula(n, cop)
  T1  <- qWeib_custom(UV[,1], xi = xi1, k = k)
  T2  <- qWeib_custom(UV[,2], xi = xi2, k = k)
  X   <- b1 * sqrt(2 * T1)
  Y   <- b2 * sqrt(2 * T2)
  cbind(X = X, Y = Y)
}

mk_rows <- function(alpha) {
  A1 <- data.frame(alpha=alpha, k=0.8, xi1=0.5, xi2=0.5, b1=0.5, b2=0.5)
  A2 <- data.frame(alpha=alpha, k=1.5, xi1=0.5, xi2=0.5, b1=0.5, b2=0.5)

  B1 <- data.frame(alpha=alpha, k=0.8, xi1=2.0, xi2=2.0, b1=2.0, b2=2.0)
  B2 <- data.frame(alpha=alpha, k=1.5, xi1=2.0, xi2=2.0, b1=2.0, b2=2.0)

  C1 <- data.frame(alpha=alpha, k=0.8, xi1=0.5, xi2=2.0, b1=0.5, b2=2.0)
  C2 <- data.frame(alpha=alpha, k=1.5, xi1=0.5, xi2=2.0, b1=0.5, b2=2.0)

  D1 <- data.frame(alpha=alpha, k=0.8, xi1=2.0, xi2=0.5, b1=2.0, b2=0.5)
  D2 <- data.frame(alpha=alpha, k=1.5, xi1=2.0, xi2=0.5, b1=2.0, b2=0.5)
  rbind(A1,A2,B1,B2,C1,C2,D1,D2)
}
params_df <- rbind(mk_rows(1.0), mk_rows(0.6), mk_rows(0.4))
row.names(params_df) <- NULL

stat_names <- c("Mean_X","Mean_Y","Mean_XY","Var_X","Var_Y","Cov_XY",
                "Skew_X","Skew_Y","Kurt_X","Kurt_Y")

results <- vector("list", nrow(params_df))
cat(sprintf("Running BWRR2 Monte Carlo (subset 24): n=%d, m=%d ...\n", n, m))

t0 <- proc.time()[3]
for (i in 1:nrow(params_df)) {
  p <- params_df[i, ]
  cat(sprintf("  [%2d/24] alpha=%.1f, k=%.1f, xi1=%.1f, xi2=%.1f, b1=%.1f, b2=%.1f\n",
              i, p$alpha, p$k, p$xi1, p$xi2, p$b1, p$b2))
  
  stats_mat <- matrix(NA_real_, nrow = m, ncol = length(stat_names))
  colnames(stats_mat) <- stat_names
  
  for (j in 1:m) {
    XY <- rBWRR2(n, alpha=p$alpha, k=p$k, xi1=p$xi1, xi2=p$xi2, b1=p$b1, b2=p$b2)
    X <- XY[,1]; Y <- XY[,2]
    stats_mat[j,"Mean_X"]  <- mean(X)
    stats_mat[j,"Mean_Y"]  <- mean(Y)
    stats_mat[j,"Mean_XY"] <- mean(X * Y)
    stats_mat[j,"Var_X"]   <- var(X)
    stats_mat[j,"Var_Y"]   <- var(Y)
    stats_mat[j,"Cov_XY"]  <- cov(X, Y)
    stats_mat[j,"Skew_X"]  <- skewness(X)
    stats_mat[j,"Skew_Y"]  <- skewness(Y)
    stats_mat[j,"Kurt_X"]  <- kurtosis(X)
    stats_mat[j,"Kurt_Y"]  <- kurtosis(Y)
  }
  
  avg_stats <- colMeans(stats_mat, na.rm = TRUE)
  results[[i]] <- cbind(params_df[i,], as.data.frame(t(avg_stats)))
}
elapsed <- proc.time()[3] - t0
cat(sprintf("Done. Elapsed: %.1f sec\n", elapsed))

out <- do.call(rbind, results)
out <- out[, c(colnames(params_df), stat_names)]

outfile <- sprintf("BWRR2_MC_subset24_n%d_m%d.csv", n, m)
write.csv(out, outfile, row.names = FALSE)
print(round(out, 4))
cat(sprintf("\nSaved: %s\n", outfile))


<pre>An Error occurred while handling another error:
yii\web\HeadersAlreadySentException: Headers already sent in  on line 0. in /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php:366
Stack trace:
#0 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php(339): yii\web\Response-&gt;sendHeaders()
#1 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/ErrorHandler.php(136): yii\web\Response-&gt;send()
#2 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/ErrorHandler.php(135): yii\web\ErrorHandler-&gt;renderException()
#3 [internal function]: yii\base\ErrorHandler-&gt;handleException()
#4 {main}
Previous exception:
yii\web\HeadersAlreadySentException: Headers already sent in  on line 0. in /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php:366
Stack trace:
#0 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Response.php(339): yii\web\Response-&gt;sendHeaders()
#1 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Application.php(656): yii\web\Response-&gt;send()
#2 /var/www/html/prof-homepages/vendor/faravaghi/yii2-filemanager/models/Files.php(696): yii\base\Application-&gt;end()
#3 /var/www/html/prof-homepages/vendor/faravaghi/yii2-filemanager/controllers/FilesController.php(484): faravaghi\filemanager\models\Files-&gt;getFile()
#4 [internal function]: faravaghi\filemanager\controllers\FilesController-&gt;actionGetFile()
#5 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/InlineAction.php(57): call_user_func_array()
#6 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Controller.php(180): yii\base\InlineAction-&gt;runWithParams()
#7 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Module.php(528): yii\base\Controller-&gt;runAction()
#8 /var/www/html/prof-homepages/vendor/yiisoft/yii2/web/Application.php(103): yii\base\Module-&gt;runAction()
#9 /var/www/html/prof-homepages/vendor/yiisoft/yii2/base/Application.php(386): yii\web\Application-&gt;handleRequest()
#10 /var/www/html/prof-homepages/backend/web/index.php(16): yii\base\Application-&gt;run()
#11 {main}</pre>