library(rgl) 

# Soucinova kopula
u <- 1:1000/1000
v <- 1:1000/1000
C=matrix(rep(0,1000*1000),nrow=1000)

for (i in 1:1000){
for (j in 1:1000){
C[i,j]<-u[i]*v[j]
}}

persp3d(u, v, C, col="skyblue",zlab='C(u,v)')
contour(u,v,C,xlab='u',ylab='v',main='Kontury C(u,v)')


# Horni kopula
for (i in 1:1000){
for (j in 1:1000){
C[i,j]<-min(u[i],v[j])
}}

persp3d(u, v, C, col="skyblue",zlab='C(u,v)')
contour(u,v,C,xlab='u',ylab='v',main='Kontury C(u,v)')


# Dolni kopula
for (i in 1:1000){
for (j in 1:1000){
C[i,j]<-max(u[i]+v[j]-1,0)
}}

persp3d(u, v, C, col="skyblue",zlab='C(u,v)')
contour(u,v,C,xlab='u',ylab='v',main='Kontury C(u,v)')


# Gumbelova kopula
theta<-2
for (i in 1:1000){
for (j in 1:1000){
C[i,j]=exp(-((-log(u[i]))^(theta)+(-log(v[j]))^(theta))^(1/theta) )
}}

persp3d(u, v, C, col="skyblue",zlab='C(u,v)')
contour(u,v,C,xlab='u',ylab='v',main='Kontury C(u,v)')


########################
library(copula)
########################
# Archimedovske kopuly #
########################

# Kopula nezavislosti
indep.cop <- indepCopula(2)

#nagenerujeme 1000 pozorovani z teto kopuly
u<-rCopula(1000, indep.cop)
plot(u,xlab='u',ylab='v')

#vykreslime kopulu a jeji kontury
par(mfrow=c(1,2))
persp (indep.cop, pCopula,xlab='u',ylab='v',main='C(u,v)')
contour(indep.cop, pCopula,nlevels=10,xlab='u',ylab='v',main='Kontury C(u,v)')

#vykreslime jeste hustotu kopuly a jeji kontury
par(mfrow=c(1,2))
persp (indep.cop, dCopula,xlab='u',ylab='v',main='c(u,v)')
contour(indep.cop, dCopula,nlevels=10,xlab='u',ylab='v',main='Kontury c(u,v)')



################################################
#Vytvareni vicerozmernych rozdeleni pomoci kopul
################################################

#pomoci kopuly a marginalnich rozdeleni definujeme vlastni distribuci

mojerozdeleni<-mvdc(gumbelCopula(3), c("unif", "gamma"),
list(list(min=0, max=100),list(shape=10, rate=1/7)))

#nagenerujeme 1000 pozorovani z tohoto rozdeleni
u<-rMvdc(1000, mojerozdeleni)
plot(u,xlab='x',ylab='y')

#vykreslime jeste distribucni funkci a sdruzenou hustotu
persp  (mojerozdeleni, pMvdc, xlim = c(0, 100), ylim=c(30, 150),xlab='x', ylab='y', main = 'F(x,y)')
persp  (mojerozdeleni, dMvdc, xlim = c(0, 100), ylim=c(30, 150),xlab='x', ylab='y', main = 'f(x,y)')
contour(mojerozdeleni, dMvdc, xlim = c(0, 100), ylim=c(30, 150),nlevels=15,xlab='x', ylab='y',main='Kontury f(x,y)')

#pridejme jeste puvodni data
points(u,col=2,pch='+')