Frank Copula

Hello all , I am trying to work on Frank Copula. Here is my artifical data set that generates correlation between x and y in R using frank copula

library(VineCopula)

Set the copula parameter
theta ← -10

Set the marginals (marginal distributions for x and y)
marginals ← c(“norm”, “norm”)

Create the Frank copula object
copula ← frankCopula(param = theta)

Generate random variables from the copula
n ← 500 # Number of data points
u ← rCopula(n, copula)

Apply inverse transformation to get simulated values for x and y
x ← qnorm(u[, 1])
y ← qnorm(u[, 2])

Plot the simulated data
plot(x, y, pch = 16, xlab = “x”, ylab = “y”, main = “Simulated Data”)

Compute and print the sample correlation
correlation ← cor(x, y)
print(paste(“Sample Correlation:”, correlation))

I am using the following Stan Code to recover the value of theta which is -10, but I am having no success, can some one please help
@spinkney I know you have done some work on this can you help please?
@andre.pfeuffer can you look into this too, please?


When including Stan code in your post it really helps if you make it as readable as possible by using Stan code chunks (```stan) with clear spacing and indentation. For example, use

functions {
real frank_copula(real u, real v, real theta) {
  return log(theta) + log1m(exp(-theta)) - theta * (u + v)
         - 2 * log(1 - exp(-theta) - (1 - exp(-theta * u)) * (1 - exp(-theta * v)));
}

}

data {
  int<lower=0> N;  // Number of observations
  vector[N] x;     // First variable
  vector[N] y;     // Second variable
}

transformed data {
  vector[N] x_scaled = (x - min(x)) / (max(x) - min(x));
  vector[N] y_scaled = (y - min(y)) / (max(y) - min(y));
}

parameters {
  real theta;  // Copula parameter
}

model {
  // Priors
  // No constraint on theta

  // Generate the copula
  vector[N] z;
  for (i in 1:N) {
    z[i] = frank_copula(x_scaled[i], y_scaled[i], theta);
  }

  // Likelihood
  for (i in 1:N) {
    target += z[i];
  }
}

generated quantities {

}
****

I believe it needs to handle the case when \theta < 0. You can easily vectorize this ldpf too (see below). I don’t think you want to do the scaling either. Just input u directly.

mod_out <- mod$sample(
  data = list(N = n,
              x = u[, 1],
              y = u[, 2]),
  parallel_chains = 4
)
functions {
real frank_copula(real u, real v, real theta) {
  return log(theta * (1 - exp(-theta))) - theta * (u + v) -
  log( (1 - exp(-theta) - (1 - exp(-theta * u)) * (1 - exp(-theta * v)))^2);
}

real frank_copula_lpdf(vector u, vector v, real theta) {
  int N = num_elements(u);
  return N * log(theta * (1 - exp(-theta))) - theta * sum(u + v) -
  sum(log( (1 - exp(-theta) - (1 - exp(-theta * u)) .* (1 - exp(-theta * v)))^2));
}

}
data {
  int<lower=0> N;  // Number of observations
  vector[N] x;     // First variable
  vector[N] y;     // Second variable
}
parameters {
  real theta;  // Copula parameter
}
model {
  // Priors
  // No constraint on theta
   target += frank_copula_lpdf(x, y, theta);
}

A post was split to a new topic: Student t copula