# Creation of a custom distribution or copulaΒΆ

In this example we are going to create a distribution or copula.

The way to go is inheriting the PythonDistribution class and overload the methods of the Distribution object.

To create a Copula, the user has to overload isCopula() and return True.

Then an instance of the new class can be passed on into a Distribution object.

At least computeCDF should be overriden.

[1]:

from __future__ import print_function
import openturns as ot
import math as m
import warnings
warnings.filterwarnings("ignore")


Inherit PythonDistribution :

[2]:

class UniformNdPy(ot.PythonDistribution):

def __init__(self, a=[0.0], b=[1.0]):
super(UniformNdPy, self).__init__(len(a))
if len(a) != len(b):
raise ValueError('Invalid bounds')
for i in range(len(a)):
if a[i] > b[i]:
raise ValueError('Invalid bounds')
self.a = a
self.b = b
self.factor = 1.0
for i in range(len(a)):
self.factor *= (b[i] - a[i])

def getRange(self):
return ot.Interval(self.a, self.b, [True] * len(self.a), [True] * len(self.a))

def getRealization(self):
X = []
for i in range(len(self.a)):
X.append(
self.a[i] + (self.b[i] - self.a[i]) * ot.RandomGenerator.Generate())
return X

def getSample(self, size):
X = []
for i in range(size):
X.append(self.getRealization())
return X

def computeCDF(self, X):
prod = 1.0
for i in range(len(self.a)):
if X[i] < self.a[i]:
return 0.0
prod *= (min(self.b[i], X[i]) - self.a[i])
return prod / self.factor

def computePDF(self, X):
for i in range(len(self.a)):
if X[i] < self.a[i]:
return 0.0
if X[i] > self.b[i]:
return 0.0
return 1.0 / self.factor

def getMean(self):
mu = []
for i in range(len(self.a)):
mu.append(0.5 * (self.a[i] + self.b[i]))
return mu

def getStandardDeviation(self):
stdev = []
for i in range(len(self.a)):
stdev.append((self.b[i] - self.a[i]) / m.sqrt(12.))
return stdev

def getSkewness(self):
return [0.] * len(self.a)

def getKurtosis(self):
return [1.8] * len(self.a)

def getStandardMoment(self, n):
if n % 2 == 1:
return [0.] * len(self.a)
return [1. / (n + 1.)] * len(self.a)

def getMoment(self, n):
return [-0.1 * n] * len(self.a)

def getCenteredMoment(self, n):
return [0.] * len(self.a)

def computeCharacteristicFunction(self, x):
if len(self.a) > 1:
raise ValueError('dim>1')
ax = self.a[0] * x
bx = self.b[0] * x
return (m.sin(bx) - m.sin(ax) + 1j * (m.cos(ax) - m.cos(bx))) / (bx - ax)

def isElliptical(self):
return (len(self.a) == 1) and (self.a[0] == -self.b[0])

def isCopula(self):
for i in range(len(self.a)):
if self.a[i] != 0.0:
return False
if self.b[i] != 1.0:
return False
return True

def getMarginal(self, indices):
subA = []
subB = []
for i in indices:
subA.append(self.a[i])
subB.append(self.b[i])
py_dist = UniformNdPy(subA, subB)
return ot.Distribution(py_dist)

def computeQuantile(self, prob, tail=False):
q = 1.0 - prob if tail else prob
quantile = self.a
for i in range(len(self.a)):
quantile[i] += q * (self.b[i] - self.a[i])
return quantile


Let us instanciate the distribution:

[3]:

distribution = ot.Distribution(UniformNdPy([5, 6], [7, 9]))


And plot the cdf:

[4]:

graph = distribution.drawCDF()
graph.setColors(["blue"])
graph

[4]:


We can easily generate sample:

[5]:

distribution.getSample(5)

[5]:

v0 v1 6.259753 8.648416 5.270553 6.097508 5.694114 8.908269 6.841359 7.50912 5.126412 6.878271

or compute the mean:

[6]:

distribution.getMean()

[6]:


[6,7.5]

Also we can compute the probability contained in an interval :

[7]:

distribution.computeProbability(ot.Interval([5.5, 6], [8.5, 9]))

[7]:

0.75


And do more (see Distribution for all methods)