# Installation des packages nécessaires (à exécuter une seule fois)
# install.packages("quantmod")

# Chargement des packages
library(quantmod)

# Fonction Black-Scholes pour évaluer les options
black_scholes <- function(S, K, r, T, sigma, type = "call") {
  # S: prix actuel de l'actif sous-jacent
  # K: prix d'exercice
  # r: taux d'intérêt sans risque (annuel)
  # T: temps jusqu'à l'échéance (en années)
  # sigma: volatilité de l'actif sous-jacent (annuelle)
  # type: "call" ou "put"
  
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  d2 <- d1 - sigma * sqrt(T)
  
  if (type == "call") {
    # Formule pour l'option d'achat (call)
    prix <- S * pnorm(d1) - K * exp(-r * T) * pnorm(d2)
  } else {
    # Formule pour l'option de vente (put)
    prix <- K * exp(-r * T) * pnorm(-d2) - S * pnorm(-d1)
  }
  
  return(prix)
}

# Définition des dates (6 mois en arrière à partir d'aujourd'hui)
date_fin <- Sys.Date()
date_debut <- date_fin - 180  # environ 6 mois

# Téléchargement des données du CAC 40
cac40 <- getSymbols("^FCHI", src = "yahoo", 
                    from = date_debut, 
                    to = date_fin, 
                    auto.assign = FALSE)

# Extraction des prix de clôture
cac40_close <- Cl(cac40)
colnames(cac40_close) <- "CAC40"

# Conversion en data frame
cac40_df <- data.frame(
  Date = index(cac40_close),
  Prix = as.numeric(cac40_close)
)

# Calcul des rendements journaliers
cac40_df <- cac40_df %>%
  arrange(Date) %>%
  mutate(Rendement = (Prix / lag(Prix) - 1) * 100)

# Calcul de la volatilité (écart-type des rendements journaliers)
volatilite_journaliere <- sd(cac40_df$Rendement, na.rm = TRUE)
volatilite_annuelle <- volatilite_journaliere * sqrt(252)  # 252 jours de trading par an

# Récupération du prix actuel (dernière valeur)
S <- tail(cac40_df$Prix, 1)
# Utilisation de la volatilité calculée
sigma <- volatilite_annuelle / 100  # Conversion en décimal

# Affichage des paramètres
cat("Prix actuel du CAC 40:", S, "\n")
cat("Volatilité annualisée:", volatilite_annuelle, "% (", sigma, ")\n")

# Paramètres pour le calcul de l'option
K <- 8100        # Prix d'exercice demandé
r <- 0.03        # Taux d'intérêt sans risque (3%, à ajuster selon le taux actuel)
T <- 3/12        # Durée jusqu'à l'échéance (3 mois = 0.25 an)

# Calcul du prix de l'option d'achat (call)
prix_call <- black_scholes(S, K, r, T, sigma, "call")

# Calcul du prix de l'option de vente (put)
prix_put <- black_scholes(S, K, r, T, sigma, "put")

# Vérification de la parité put-call
parite <- prix_call + K * exp(-r * T) - prix_put - S

# Calcul des Grecques
# Delta (sensibilité du prix de l'option par rapport au prix de l'actif sous-jacent)
delta_call <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  return(pnorm(d1))
}

# Delta put
delta_put <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  return(pnorm(d1) - 1)
}

# Gamma (taux de variation du delta)
gamma <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  return(dnorm(d1) / (S * sigma * sqrt(T)))
}

# Vega (sensibilité à la volatilité)
vega <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  return(S * dnorm(d1) * sqrt(T) / 100)  # divisé par 100 pour exprimer en % de volatilité
}

# Thêta (sensibilité au temps)
theta_call <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  d2 <- d1 - sigma * sqrt(T)
  return(-S * dnorm(d1) * sigma / (2 * sqrt(T)) - r * K * exp(-r * T) * pnorm(d2))
}

# Thêta put
theta_put <- function(S, K, r, T, sigma) {
  d1 <- (log(S/K) + (r + sigma^2/2) * T) / (sigma * sqrt(T))
  d2 <- d1 - sigma * sqrt(T)
  return(-S * dnorm(d1) * sigma / (2 * sqrt(T)) + r * K * exp(-r * T) * pnorm(-d2))
}

# Calcul des grecques pour nos options
delta_c <- delta_call(S, K, r, T, sigma)
delta_p <- delta_put(S, K, r, T, sigma)
gamma_val <- gamma(S, K, r, T, sigma)  # Même gamma pour call et put
vega_val <- vega(S, K, r, T, sigma)    # Même vega pour call et put
theta_c <- theta_call(S, K, r, T, sigma) / 365  # Convertir en valeur journalière
theta_p <- theta_put(S, K, r, T, sigma) / 365   # Convertir en valeur journalière

# Affichage des résultats
cat("\n----- RÉSULTATS DES OPTIONS SUR LE CAC 40 -----\n")
cat("Prix d'exercice (Strike):", K, "\n")
cat("Temps jusqu'à l'échéance:", T*12, "mois\n")
cat("Taux sans risque:", r*100, "%\n\n")
cat("Prix de l'option d'achat (call):", round(prix_call, 2), "EUR\n")
cat("Prix de l'option de vente (put):", round(prix_put, 2), "EUR\n")
cat("Vérification parité put-call (doit être proche de 0):", round(parite, 4), "\n\n")
cat("------- GRECQUES -------\n")
cat("Delta call:", round(delta_c, 4), "\n")
cat("Delta put:", round(delta_p, 4), "\n")
cat("Gamma (identique call/put):", round(gamma_val, 6), "\n")
cat("Vega (identique call/put):", round(vega_val, 2), "(pour 1% de volatilité)\n")
cat("Theta call:", round(theta_c, 2), "par jour\n")
cat("Theta put:", round(theta_p, 2), "par jour\n")