1 // -*- C++ -*- 2 /** 3 * @brief Abstract top-level class for all distributions 4 * 5 * Copyright 2005-2021 Airbus-EDF-IMACS-ONERA-Phimeca 6 * 7 * This library is free software: you can redistribute it and/or modify 8 * it under the terms of the GNU Lesser General Public License as published by 9 * the Free Software Foundation, either version 3 of the License, or 10 * (at your option) any later version. 11 * 12 * This library is distributed in the hope that it will be useful, 13 * but WITHOUT ANY WARRANTY; without even the implied warranty of 14 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE. See the 15 * GNU Lesser General Public License for more details. 16 * 17 * You should have received a copy of the GNU Lesser General Public License 18 * along with this library. If not, see <http://www.gnu.org/licenses/>. 19 * 20 */ 21 #ifndef OPENTURNS_DISTRIBUTIONIMPLEMENTATION_HXX 22 #define OPENTURNS_DISTRIBUTIONIMPLEMENTATION_HXX 23 24 #include "openturns/PersistentObject.hxx" 25 #include "openturns/Point.hxx" 26 #include "openturns/PointWithDescription.hxx" 27 #include "openturns/Sample.hxx" 28 #include "openturns/Indices.hxx" 29 #include "openturns/Interval.hxx" 30 #include "openturns/LevelSet.hxx" 31 #include "openturns/CorrelationMatrix.hxx" 32 #include "openturns/SquareMatrix.hxx" 33 #include "openturns/Graph.hxx" 34 #include "openturns/Description.hxx" 35 #include "openturns/EvaluationImplementation.hxx" 36 #include "openturns/GradientImplementation.hxx" 37 #include "openturns/UniVariateFunctionImplementation.hxx" 38 #include "openturns/PersistentCollection.hxx" 39 #include "openturns/UniVariatePolynomial.hxx" 40 #include "openturns/PiecewiseHermiteEvaluation.hxx" 41 #include "openturns/ResourceMap.hxx" 42 43 BEGIN_NAMESPACE_OPENTURNS 44 45 // Forward declaration 46 class Distribution; 47 48 /** 49 * @class DistributionImplementation 50 * 51 * The class describes the probabilistic concept of distribution. 52 * Instances of this class can produce samples following the 53 * distribution, can compute PDF or CDF, etc. 54 * They are the actual key component of RandomVectors. 55 */ 56 class OT_API DistributionImplementation 57 : public PersistentObject 58 { 59 CLASSNAME 60 public: 61 62 typedef Pointer<DistributionImplementation> Implementation; 63 typedef Function IsoProbabilisticTransformation; 64 typedef IsoProbabilisticTransformation InverseIsoProbabilisticTransformation; 65 typedef Collection<Point> PointCollection; 66 typedef Collection<PointWithDescription> PointWithDescriptionCollection; 67 68 /** Default constructor */ 69 DistributionImplementation(); 70 71 /** Comparison operator */ 72 Bool operator ==(const DistributionImplementation & other) const; 73 protected: 74 virtual Bool equals(const DistributionImplementation & other) const; 75 public: 76 Bool operator !=(const DistributionImplementation & other) const; 77 78 /** Addition operator */ 79 Distribution operator + (const DistributionImplementation & other) const; 80 Distribution operator + (const Distribution & other) const; 81 82 Distribution operator + (const Scalar value) const; 83 84 /** Subtraction operator */ 85 Distribution operator - (const DistributionImplementation & other) const; 86 Distribution operator - (const Distribution & other) const; 87 88 Distribution operator - (const Scalar value) const; 89 90 /** Multiplication operator */ 91 Distribution operator * (const DistributionImplementation & other) const; 92 Distribution operator * (const Distribution & other) const; 93 94 Distribution operator * (const Scalar value) const; 95 96 /** Division operator */ 97 Distribution operator / (const DistributionImplementation & other) const; 98 Distribution operator / (const Distribution & other) const; 99 100 Distribution operator / (const Scalar value) const; 101 102 /** Methods to transform distributions by usual functions */ 103 Distribution cos() const; 104 Distribution sin() const; 105 Distribution tan() const; 106 107 Distribution acos() const; 108 Distribution asin() const; 109 Distribution atan() const; 110 111 Distribution cosh() const; 112 Distribution sinh() const; 113 Distribution tanh() const; 114 115 Distribution acosh() const; 116 Distribution asinh() const; 117 Distribution atanh() const; 118 119 Distribution exp() const; 120 Distribution log() const; 121 Distribution ln() const; 122 123 Distribution pow(const SignedInteger exponent) const; 124 Distribution pow(const Scalar exponent) const; 125 Distribution inverse() const; 126 Distribution sqr() const; 127 Distribution sqrt() const; 128 Distribution cbrt() const; 129 Distribution abs() const; 130 131 /** String converter */ 132 String __repr__() const override; 133 String __str__(const String & offset = "") const override; 134 135 /** Weight accessor */ 136 void setWeight(Scalar w); 137 Scalar getWeight() const; 138 139 /** Dimension accessor */ 140 UnsignedInteger getDimension() const; 141 142 /* Here is the interface that all derived class must implement */ 143 144 /** Virtual constructor */ 145 DistributionImplementation * clone() const override; 146 147 /** Get one realization of the distribution */ 148 virtual Point getRealization() const; 149 protected: 150 virtual Point getRealizationByInversion() const; 151 152 public: 153 /** Get a numerical sample whose elements follow the distributionImplementation */ 154 virtual Sample getSample(const UnsignedInteger size) const; 155 virtual Sample getSampleByInversion(const UnsignedInteger size) const; 156 virtual Sample getSampleByQMC(const UnsignedInteger size) const; 157 158 /** Get the DDF of the distribution */ 159 virtual Scalar computeDDF(const Scalar scalar) const; 160 virtual Point computeDDF(const Point & point) const; 161 virtual Sample computeDDF(const Sample & sample) const; 162 protected: 163 virtual Sample computeDDFSequential(const Sample & sample) const; 164 virtual Sample computeDDFParallel(const Sample & sample) const; 165 public: 166 167 /** Get the PDF of the distribution */ 168 virtual Scalar computePDF(const Scalar scalar) const; 169 virtual Scalar computePDF(const Point & point) const; 170 virtual Sample computePDF(const Sample & sample) const; 171 protected: 172 virtual Sample computePDFSequential(const Sample & sample) const; 173 virtual Sample computePDFParallel(const Sample & sample) const; 174 public: 175 176 virtual Scalar computeLogPDF(const Scalar scalar) const; 177 virtual Scalar computeLogPDF(const Point & point) const; 178 virtual Sample computeLogPDF(const Sample & sample) const; 179 protected: 180 virtual Sample computeLogPDFSequential(const Sample & sample) const; 181 virtual Sample computeLogPDFParallel(const Sample & sample) const; 182 public: 183 184 /** Compute the PDF of 1D distributions over a regular grid */ 185 virtual Sample computePDF(const Scalar xMin, 186 const Scalar xMax, 187 const UnsignedInteger pointNumber, 188 Sample & gridOut) const; 189 190 /** Compute the PDF of nD distributions over a regular grid */ 191 virtual Sample computePDF(const Point & xMin, 192 const Point & xMax, 193 const Indices & pointNumber, 194 Sample & gridOut) const; 195 196 /** Compute the log-PDF of 1D distributions over a regular grid */ 197 virtual Sample computeLogPDF(const Scalar xMin, 198 const Scalar xMax, 199 const UnsignedInteger pointNumber, 200 Sample & gridOut) const; 201 202 /** Compute the log-PDF of nD distributions over a regular grid */ 203 virtual Sample computeLogPDF(const Point & xMin, 204 const Point & xMax, 205 const Indices & pointNumber, 206 Sample & gridOut) const; 207 208 /** Get the CDF of the distribution */ 209 virtual Scalar computeCDF(const Scalar scalar) const; 210 virtual Scalar computeComplementaryCDF(const Scalar scalar) const; 211 virtual Scalar computeSurvivalFunction(const Scalar scalar) const; 212 213 virtual Scalar computeCDF(const Point & point) const; 214 virtual Scalar computeComplementaryCDF(const Point & point) const; 215 virtual Scalar computeSurvivalFunction(const Point & point) const; 216 virtual Point computeInverseSurvivalFunction(const Scalar point) const; 217 #ifndef SWIG 218 virtual Point computeInverseSurvivalFunction(const Scalar prob, 219 Scalar & marginalProb) const; 220 #endif 221 protected: 222 virtual Sample computeCDFSequential(const Sample & sample) const; 223 virtual Sample computeCDFParallel(const Sample & sample) const; 224 public: 225 virtual Sample computeCDF(const Sample & sample) const; 226 protected: 227 virtual Sample computeSurvivalFunctionSequential(const Sample & sample) const; 228 virtual Sample computeSurvivalFunctionParallel(const Sample & sample) const; 229 public: 230 virtual Sample computeSurvivalFunction(const Sample & sample) const; 231 protected: 232 virtual Sample computeComplementaryCDFSequential(const Sample & sample) const; 233 virtual Sample computeComplementaryCDFParallel(const Sample & sample) const; 234 public: 235 virtual Sample computeComplementaryCDF(const Sample & sample) const; 236 237 /** Compute the CDF of 1D distributions over a regular grid */ 238 virtual Sample computeCDF(const Scalar xMin, 239 const Scalar xMax, 240 const UnsignedInteger pointNumber, 241 Sample & gridOut) const; 242 243 /** Compute the CDF of nD distributions over a regular grid */ 244 virtual Sample computeCDF(const Point & xMin, 245 const Point & xMax, 246 const Indices & pointNumber, 247 Sample & gridOut) const; 248 249 virtual Sample computeComplementaryCDF(const Scalar xMin, 250 const Scalar xMax, 251 const UnsignedInteger pointNumber, 252 Sample & gridOut) const; 253 254 /** Get the probability content of an interval */ 255 virtual Scalar computeProbability(const Interval & interval) const; 256 protected: 257 /** Generic implementation for continuous distributions */ 258 virtual Scalar computeProbabilityContinuous(const Interval & interval) const; 259 /** Generic implementation for discrete distributions */ 260 virtual Scalar computeProbabilityDiscrete(const Interval & interval) const; 261 /** Generic implementation for general distributions */ 262 virtual Scalar computeProbabilityGeneral(const Interval & interval) const; 263 public: 264 265 /** Get the characteristic function of the distribution, i.e. phi(u) = E(exp(I*u*X)) */ 266 virtual Complex computeCharacteristicFunction(const Scalar x) const; 267 virtual Complex computeCharacteristicFunction(const Point & x) const; 268 virtual Complex computeLogCharacteristicFunction(const Scalar x) const; 269 virtual Complex computeLogCharacteristicFunction(const Point & x) const; 270 271 /** Compute the generating function, i.e. psi(z) = E(z^X) */ 272 virtual Scalar computeGeneratingFunction(const Scalar z) const; 273 virtual Scalar computeLogGeneratingFunction(const Scalar z) const; 274 275 virtual Complex computeGeneratingFunction(const Complex & z) const; 276 virtual Complex computeLogGeneratingFunction(const Complex & z) const; 277 278 /** Compute the entropy of the distribution */ 279 virtual Scalar computeEntropy() const; 280 281 /** Get the PDF gradient of the distribution */ 282 virtual Point computePDFGradient(const Point & point) const; 283 virtual Sample computePDFGradient(const Sample & inSample) const; 284 285 /** Get the logPDF gradient of the distribution */ 286 virtual Point computeLogPDFGradient(const Point & point) const; 287 virtual Sample computeLogPDFGradient(const Sample & inSample) const; 288 289 protected: 290 virtual Sample computeLogPDFGradientSequential(const Sample & sample) const; 291 virtual Sample computeLogPDFGradientParallel(const Sample & sample) const; 292 293 public: 294 /** Get the CDF gradient of the distribution */ 295 virtual Point computeCDFGradient(const Point & point) const; 296 virtual Sample computeCDFGradient(const Sample & inSample) const; 297 public: 298 299 /** Get the quantile of the distribution */ 300 virtual Point computeQuantile(const Scalar prob, 301 const Bool tail = false) const; 302 #ifndef SWIG 303 virtual Point computeQuantile(const Scalar prob, 304 const Bool tail, 305 Scalar & marginalProb) const; 306 #endif 307 /** Quantile computation for dimension=1 */ 308 virtual Scalar computeScalarQuantile(const Scalar prob, 309 const Bool tail = false) const; 310 311 /** Get the quantile over a provided grid */ 312 protected: 313 virtual Sample computeQuantileSequential(const Point & prob, 314 const Bool tail = false) const; 315 virtual Sample computeQuantileParallel(const Point & prob, 316 const Bool tail = false) const; 317 Point computeQuantileCopula(const Scalar prob, const Bool tail = false) const; 318 public: 319 virtual Sample computeQuantile(const Point & prob, 320 const Bool tail = false) const; 321 322 /** Compute the quantile over a regular grid */ 323 virtual Sample computeQuantile(const Scalar qMin, 324 const Scalar qMax, 325 const UnsignedInteger pointNumber, 326 const Bool tail = false) const; 327 328 #ifndef SWIG 329 virtual Sample computeQuantile(const Scalar qMin, 330 const Scalar qMax, 331 const UnsignedInteger pointNumber, 332 Sample & grid, 333 const Bool tail = false) const; 334 #endif 335 336 /** Get the product minimum volume interval containing a given probability of the distribution */ 337 virtual Interval computeMinimumVolumeInterval(const Scalar prob) const; 338 virtual Interval computeMinimumVolumeIntervalWithMarginalProbability(const Scalar prob, Scalar & marginalProbOut) const; 339 340 protected: 341 Interval computeUnivariateMinimumVolumeIntervalByOptimization(const Scalar prob, 342 Scalar & marginalProbOut) const; 343 Interval computeUnivariateMinimumVolumeIntervalByRootFinding(const Scalar prob, 344 Scalar & marginalProbOut) const; 345 346 public: 347 348 /** Get the product bilateral confidence interval containing a given probability of the distribution */ 349 virtual Interval computeBilateralConfidenceInterval(const Scalar prob) const; 350 virtual Interval computeBilateralConfidenceIntervalWithMarginalProbability(const Scalar prob, Scalar & marginalProbOut) const; 351 352 /** Get the product unilateral confidence interval containing a given probability of the distribution */ 353 virtual Interval computeUnilateralConfidenceInterval(const Scalar prob, const Bool tail = false) const; 354 virtual Interval computeUnilateralConfidenceIntervalWithMarginalProbability(const Scalar prob, const Bool tail, Scalar & marginalProbOut) const; 355 356 /** Get the minimum volume level set containing a given probability of the distribution */ 357 virtual LevelSet computeMinimumVolumeLevelSet(const Scalar prob) const; 358 virtual LevelSet computeMinimumVolumeLevelSetWithThreshold(const Scalar prob, Scalar & thresholdOut) const; 359 360 protected: 361 virtual LevelSet computeUnivariateMinimumVolumeLevelSetByQMC(const Scalar prob, 362 Scalar & thresholdOut) const; 363 364 public: 365 366 /** Get the mathematical and numerical range of the distribution. 367 Its mathematical range is the smallest closed interval outside 368 of which the PDF is zero, and the numerical range is the interval 369 outside of which the PDF is rounded to zero in double precision */ 370 virtual Interval getRange() const; 371 372 protected: 373 374 virtual void setRange(const Interval & range); 375 376 public: 377 /** Get the roughness, i.e. the L2-norm of the PDF */ 378 virtual Scalar getRoughness() const; 379 380 /** Get the mean of the distribution */ 381 virtual Point getMean() const; 382 383 /** Get the standard deviation of the distribution */ 384 virtual Point getStandardDeviation() const; 385 386 /** Get the skewness of the distribution */ 387 virtual Point getSkewness() const; 388 389 /** Get the kurtosis of the distribution */ 390 virtual Point getKurtosis() const; 391 392 /** Get the raw moments of the standardized distribution */ 393 virtual Point getStandardMoment(const UnsignedInteger n) const; 394 395 /** Get the raw moments of the distribution */ 396 virtual Point getMoment(const UnsignedInteger n) const; 397 398 /** Get the centered moments of the distribution */ 399 virtual Point getCenteredMoment(const UnsignedInteger n) const; 400 401 /** Get the shifted moments of the distribution */ 402 virtual Point getShiftedMoment(const UnsignedInteger n, 403 const Point & shift) const; 404 405 /** Get the covariance of the distribution */ 406 virtual CovarianceMatrix getCovariance() const; 407 408 /** Correlation matrix accessor */ 409 CorrelationMatrix getCorrelation() const; 410 411 /** Get the Pearson correlation of the distribution */ 412 virtual CorrelationMatrix getPearsonCorrelation() const; 413 414 /** Get the Spearman correlation of the distribution */ 415 virtual CorrelationMatrix getSpearmanCorrelation() const; 416 417 /** Get the Kendall concordance of the distribution */ 418 virtual CorrelationMatrix getKendallTau() const; 419 420 /** Get the shape matrix of the distribution, ie the correlation matrix 421 of its copula if it is elliptical */ 422 virtual CorrelationMatrix getShapeMatrix() const; 423 424 /** Cholesky factor of the covariance matrix accessor */ 425 TriangularMatrix getCholesky() const; 426 427 /** Inverse of the Cholesky factor of the covariance matrix accessor */ 428 TriangularMatrix getInverseCholesky() const; 429 430 /** Check if the distribution is a copula */ 431 virtual Bool isCopula() const; 432 433 /** Check if the distribution is elliptical */ 434 virtual Bool isElliptical() const; 435 436 /** Check if the distribution is continuous */ 437 virtual Bool isContinuous() const; 438 439 /** Check if the distribution is discrete */ 440 virtual Bool isDiscrete() const; 441 442 /** Tell if the distribution is integer valued */ 443 virtual Bool isIntegral() const; 444 445 /** Tell if the distribution has elliptical copula */ 446 virtual Bool hasEllipticalCopula() const; 447 448 /** Tell if the distribution has independent copula */ 449 virtual Bool hasIndependentCopula() const; 450 451 /** Get the support of a distribution that intersect a given interval */ 452 virtual Sample getSupport(const Interval & interval) const; 453 454 /** Get the support on the whole range */ 455 virtual Sample getSupport() const; 456 457 /** Get the discrete probability levels */ 458 virtual Point getProbabilities() const; 459 460 /** Get the PDF singularities inside of the range - 1D only */ 461 virtual Point getSingularities() const; 462 463 /** Compute the density generator of the elliptical generator, i.e. 464 * the function phi such that the density of the distribution can 465 * be written as p(x) = phi(t(x-mu)R(x-mu)) */ 466 virtual Scalar computeDensityGenerator(const Scalar betaSquare) const; 467 468 /** Compute the derivative of the density generator */ 469 virtual Scalar computeDensityGeneratorDerivative(const Scalar betaSquare) const; 470 471 /** Compute the seconde derivative of the density generator */ 472 virtual Scalar computeDensityGeneratorSecondDerivative(const Scalar betaSquare) const; 473 474 /** Compute the radial distribution CDF */ 475 virtual Scalar computeRadialDistributionCDF(const Scalar radius, 476 const Bool tail = false) const; 477 478 /** Get the i-th marginal distribution */ 479 virtual Distribution getMarginal(const UnsignedInteger i) const; 480 481 /** Get the distribution of the marginal distribution corresponding to indices dimensions */ 482 virtual Distribution getMarginal(const Indices & indices) const; 483 484 /** Get the copula of a distribution */ 485 virtual Distribution getCopula() const; 486 487 /** Compute the DDF of Xi | X1, ..., Xi-1. x = Xi, y = (X1,...,Xi-1) */ 488 virtual Scalar computeConditionalDDF(const Scalar x, 489 const Point & y) const; 490 491 virtual Point computeSequentialConditionalDDF(const Point & x) const; 492 493 /** Compute the PDF of Xi | X1, ..., Xi-1. x = Xi, y = (X1,...,Xi-1) */ 494 virtual Scalar computeConditionalPDF(const Scalar x, 495 const Point & y) const; 496 497 virtual Point computeSequentialConditionalPDF(const Point & x) const; 498 499 virtual Point computeConditionalPDF(const Point & x, 500 const Sample & y) const; 501 502 /** Compute the CDF of Xi | X1, ..., Xi-1. x = Xi, y = (X1,...,Xi-1) */ 503 virtual Scalar computeConditionalCDF(const Scalar x, 504 const Point & y) const; 505 506 virtual Point computeSequentialConditionalCDF(const Point & x) const; 507 508 virtual Point computeConditionalCDF(const Point & x, 509 const Sample & y) const; 510 511 /** Compute the quantile of Xi | X1, ..., Xi-1, i.e. x such that CDF(x|y) = q with x = Xi, y = (X1,...,Xi-1) */ 512 virtual Scalar computeConditionalQuantile(const Scalar q, 513 const Point & y) const; 514 515 virtual Point computeSequentialConditionalQuantile(const Point & q) const; 516 517 virtual Point computeConditionalQuantile(const Point & q, 518 const Sample & y) const; 519 520 /** Get the isoprobabilist transformation */ 521 virtual IsoProbabilisticTransformation getIsoProbabilisticTransformation() const; 522 523 /** Get the inverse isoprobabilist transformation */ 524 virtual InverseIsoProbabilisticTransformation getInverseIsoProbabilisticTransformation() const; 525 526 /** Get the standard distribution */ 527 virtual Distribution getStandardDistribution() const; 528 529 /** Get the standard representative in the parametric family, associated with the standard moments */ 530 virtual Distribution getStandardRepresentative() const; 531 532 /** integrationNodesNumber accessors */ 533 UnsignedInteger getIntegrationNodesNumber() const; 534 void setIntegrationNodesNumber(const UnsignedInteger integrationNodesNumber) const; 535 536 protected: 537 /** Gauss nodes and weights accessor */ 538 Point getGaussNodesAndWeights(Point & weights) const; 539 540 public: 541 /** Draw the PDF of the distribution when its dimension is 1 or 2 */ 542 virtual Graph drawPDF(const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 543 const Bool logScale = false) const; 544 545 /** Draw the PDF of the distribution when its dimension is 1 */ 546 virtual Graph drawPDF(const Scalar xMin, 547 const Scalar xMax, 548 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 549 const Bool logScale = false) const; 550 551 /** Draw the PDF of a 1D marginal */ 552 virtual Graph drawMarginal1DPDF(const UnsignedInteger marginalIndex, 553 const Scalar xMin, 554 const Scalar xMax, 555 const UnsignedInteger pointNumber, 556 const Bool logScale = false) const; 557 558 /** Draw the PDF of the distribution when its dimension is 2 */ 559 virtual Graph drawPDF(const Point & xMin, 560 const Point & xMax, 561 const Indices & pointNumber, 562 const Bool logScaleX = false, 563 const Bool logScaleY = false) const; 564 virtual Graph drawPDF(const Point & xMin, 565 const Point & xMax, 566 const Bool logScaleX = false, 567 const Bool logScaleY = false) const; 568 virtual Graph drawPDF(const Indices & pointNumber, 569 const Bool logScaleX = false, 570 const Bool logScaleY = false) const; 571 572 /** Draw the PDF of a 2D marginal */ 573 virtual Graph drawMarginal2DPDF(const UnsignedInteger firstMarginal, 574 const UnsignedInteger secondMarginal, 575 const Point & xMin, 576 const Point & xMax, 577 const Indices & pointNumber, 578 const Bool logScaleX = false, 579 const Bool logScaleY = false) const; 580 581 /** Draw the log-PDF of the distribution when its dimension is 1 or 2 */ 582 virtual Graph drawLogPDF(const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 583 const Bool logScale = false) const; 584 585 /** Draw the log-PDF of the distribution when its dimension is 1 */ 586 virtual Graph drawLogPDF(const Scalar xMin, 587 const Scalar xMax, 588 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 589 const Bool logScale = false) const; 590 591 /** Draw the log-PDF of a 1D marginal */ 592 virtual Graph drawMarginal1DLogPDF(const UnsignedInteger marginalIndex, 593 const Scalar xMin, 594 const Scalar xMax, 595 const UnsignedInteger pointNumber, 596 const Bool logScale = false) const; 597 598 /** Draw the log-PDF of the distribution when its dimension is 2 */ 599 virtual Graph drawLogPDF(const Point & xMin, 600 const Point & xMax, 601 const Indices & pointNumber, 602 const Bool logScaleX = false, 603 const Bool logScaleY = false) const; 604 virtual Graph drawLogPDF(const Point & xMin, 605 const Point & xMax, 606 const Bool logScaleX = false, 607 const Bool logScaleY = false) const; 608 virtual Graph drawLogPDF(const Indices & pointNumber, 609 const Bool logScaleX = false, 610 const Bool logScaleY = false) const; 611 612 /** Draw the PDF of a 2D marginal */ 613 virtual Graph drawMarginal2DLogPDF(const UnsignedInteger firstMarginal, 614 const UnsignedInteger secondMarginal, 615 const Point & xMin, 616 const Point & xMax, 617 const Indices & pointNumber, 618 const Bool logScaleX = false, 619 const Bool logScaleY = false) const; 620 621 /** Draw the CDF of the distribution when its dimension is 1 or 2 */ 622 virtual Graph drawCDF(const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 623 const Bool logScale = false) const; 624 625 /** Draw the CDF of the distribution when its dimension is 1 */ 626 virtual Graph drawCDF(const Scalar xMin, 627 const Scalar xMax, 628 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 629 const Bool logScale = false) const; 630 631 /** Draw the CDF of the distribution when its dimension is 2 */ 632 virtual Graph drawCDF(const Point & xMin, 633 const Point & xMax, 634 const Indices & pointNumber, 635 const Bool logScaleX = false, 636 const Bool logScaleY = false) const; 637 virtual Graph drawCDF(const Point & xMin, 638 const Point & xMax, 639 const Bool logScaleX = false, 640 const Bool logScaleY = false) const; 641 virtual Graph drawCDF(const Indices & pointNumber, 642 const Bool logScaleX = false, 643 const Bool logScaleY = false) const; 644 645 /** Draw the CDF of a 1D marginal */ 646 virtual Graph drawMarginal1DCDF(const UnsignedInteger marginalIndex, 647 const Scalar xMin, 648 const Scalar xMax, 649 const UnsignedInteger pointNumber, 650 const Bool logScale = false) const; 651 652 /** Draw the CDF of a 2D marginal */ 653 virtual Graph drawMarginal2DCDF(const UnsignedInteger firstMarginal, 654 const UnsignedInteger secondMarginal, 655 const Point & xMin, 656 const Point & xMax, 657 const Indices & pointNumber, 658 const Bool logScaleX = false, 659 const Bool logScaleY = false) const; 660 661 /** Draw the SurvivalFunction of the distribution when its dimension is 1 or 2 */ 662 virtual Graph drawSurvivalFunction(const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 663 const Bool logScale = false) const; 664 665 /** Draw the SurvivalFunction of the distribution when its dimension is 1 */ 666 virtual Graph drawSurvivalFunction(const Scalar xMin, 667 const Scalar xMax, 668 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 669 const Bool logScale = false) const; 670 671 /** Draw the SurvivalFunction of the distribution when its dimension is 2 */ 672 virtual Graph drawSurvivalFunction(const Point & xMin, 673 const Point & xMax, 674 const Indices & pointNumber, 675 const Bool logScaleX = false, 676 const Bool logScaleY = false) const; 677 virtual Graph drawSurvivalFunction(const Point & xMin, 678 const Point & xMax, 679 const Bool logScaleX = false, 680 const Bool logScaleY = false) const; 681 virtual Graph drawSurvivalFunction(const Indices & pointNumber, 682 const Bool logScaleX = false, 683 const Bool logScaleY = false) const; 684 685 /** Draw the SurvivalFunction of a 1D marginal */ 686 virtual Graph drawMarginal1DSurvivalFunction(const UnsignedInteger marginalIndex, 687 const Scalar xMin, 688 const Scalar xMax, 689 const UnsignedInteger pointNumber, 690 const Bool logScale = false) const; 691 692 /** Draw the SurvivalFunction of a 2D marginal */ 693 virtual Graph drawMarginal2DSurvivalFunction(const UnsignedInteger firstMarginal, 694 const UnsignedInteger secondMarginal, 695 const Point & xMin, 696 const Point & xMax, 697 const Indices & pointNumber, 698 const Bool logScaleX = false, 699 const Bool logScaleY = false) const; 700 701 /** Draw the quantile of the distribution when its dimension is 1 or 2 */ 702 virtual Graph drawQuantile(const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 703 const Bool logScale = false) const; 704 705 virtual Graph drawQuantile(const Scalar qMin, 706 const Scalar qMax, 707 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 708 const Bool logScale = false) const; 709 710 protected: 711 virtual Graph drawQuantile1D(const Scalar qMin, 712 const Scalar qMax, 713 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 714 const Bool logScale = false) const; 715 716 virtual Graph drawQuantile2D(const Scalar qMin, 717 const Scalar qMax, 718 const UnsignedInteger pointNumber = ResourceMap::GetAsUnsignedInteger("Distribution-DefaultPointNumber"), 719 const Bool logScaleX = false, 720 const Bool logScaleY = false) const; 721 public: 722 723 /** Parameters value and description accessor */ 724 virtual PointWithDescriptionCollection getParametersCollection() const; 725 virtual void setParametersCollection(const PointWithDescriptionCollection & parametersCollection); 726 virtual void setParametersCollection(const PointCollection & parametersCollection); 727 728 /** Parameters value accessor */ 729 virtual Point getParameter() const; 730 virtual void setParameter(const Point & parameters); 731 732 /** Parameters description accessor */ 733 virtual Description getParameterDescription() const; 734 735 /** Total number of parameters */ 736 virtual UnsignedInteger getParameterDimension() const; 737 738 /** Description accessor */ 739 void setDescription(const Description & description); 740 Description getDescription() const; 741 742 /** Method save() stores the object through the StorageManager */ 743 void save(Advocate & adv) const override; 744 745 /** Method load() reloads the object from the StorageManager */ 746 void load(Advocate & adv) override; 747 748 /** Accessor to PDF computation precision */ 749 Scalar getPDFEpsilon() const; 750 751 /** Accessor to CDF computation precision */ 752 Scalar getCDFEpsilon() const; 753 754 /** Get a positon indicator for a 1D distribution */ 755 Scalar getPositionIndicator() const; 756 757 /** Get a dispersion indicator for a 1D distribution */ 758 Scalar getDispersionIndicator() const; 759 760 /** Is it safe to compute PDF/CDF etc in parallel? */ 761 Bool isParallel() const; 762 void setParallel(const Bool flag); 763 764 /** Quantile epsilon accessor */ 765 Scalar getQuantileEpsilon() const; 766 void setQuantileEpsilon(const Scalar quantileEpsilon); 767 768 protected: 769 770 /** Draw the PDF of a discrete distribution */ 771 virtual Graph drawDiscretePDF(const Scalar xMin, 772 const Scalar xMax, 773 const Bool logScale = false) const; 774 775 /** Draw the log-PDF of a discrete distribution */ 776 virtual Graph drawDiscreteLogPDF(const Scalar xMin, 777 const Scalar xMax, 778 const Bool logScale = false) const; 779 780 /** Draw the CDF of a discrete distribution */ 781 Graph drawDiscreteCDF(const Scalar xMin, 782 const Scalar xMax, 783 const Bool logScale = false) const; 784 785 /** Draw the SurvivalFunction of a discrete distribution */ 786 Graph drawDiscreteSurvivalFunction(const Scalar xMin, 787 const Scalar xMax, 788 const Bool logScale = false) const; 789 790 /** Compute the characteristic function of 1D distributions in a regular pattern with cache */ 791 virtual Complex computeCharacteristicFunction(const UnsignedInteger index, 792 const Scalar step) const; 793 virtual Complex computeLogCharacteristicFunction(const UnsignedInteger index, 794 const Scalar step) const; 795 virtual Complex computeCharacteristicFunction(const Indices & indices, 796 const Point & step) const; 797 virtual Complex computeLogCharacteristicFunction(const Indices & indices, 798 const Point & step) const; 799 800 /** Compute the mean of the distribution */ 801 virtual void computeMean() const; 802 803 /** Compute the covariance of the distribution */ 804 virtual void computeCovariance() const; 805 virtual void computeCovarianceCopula() const; 806 virtual void computeCovarianceContinuous() const; 807 virtual void computeCovarianceDiscrete() const; 808 virtual void computeCovarianceGeneral() const; 809 810 /** Compute the shifted moments of the distribution */ 811 virtual Point computeShiftedMomentContinuous(const UnsignedInteger n, 812 const Point & shift) const; 813 virtual Point computeShiftedMomentDiscrete(const UnsignedInteger n, 814 const Point & shift) const; 815 virtual Point computeShiftedMomentGeneral(const UnsignedInteger n, 816 const Point & shift) const; 817 818 819 /** Compute the nodes and weights of the 1D gauss integration rule over [-1, 1] */ 820 virtual void computeGaussNodesAndWeights() const; 821 822 /** Dimension accessor */ 823 void setDimension(UnsignedInteger dim); 824 825 /** Compute the numerical range of the distribution given the parameters values */ 826 virtual void computeRange(); 827 virtual Point computeLowerBound() const; 828 virtual Point computeUpperBound() const; 829 830 /** Compute the standard distribution associated with the current distribution */ 831 virtual void computeStandardDistribution() const; 832 833 /** Interpolate the PDF and CDF for smooth continuous distributions */ 834 virtual Collection<PiecewiseHermiteEvaluation> interpolatePDFCDF(const UnsignedInteger n); 835 836 mutable Point mean_; 837 mutable CovarianceMatrix covariance_; 838 mutable Point gaussNodes_; 839 mutable Point gaussWeights_; 840 841 /** The integration nodes number for covariance computation */ 842 mutable UnsignedInteger integrationNodesNumber_; 843 844 /** Flags to avoid redundant computations */ 845 mutable Bool isAlreadyComputedMean_; 846 mutable Bool isAlreadyComputedCovariance_; 847 mutable Bool isAlreadyComputedGaussNodesAndWeights_; 848 849 /** Indicators for PDF and CDF absolute precision. By default, equals to DefaultPDFEpsilon, DefaultCDFEpsilon and DefaultQuantileEpsilon */ 850 Scalar pdfEpsilon_; 851 Scalar cdfEpsilon_; 852 Scalar quantileEpsilon_; 853 UnsignedInteger quantileIterations_; 854 855 /** Standard distribution */ 856 mutable Bool isAlreadyComputedStandardDistribution_; 857 // We cannot have a Distribution member in DistributionImplementation, there would 858 // be circular dependency. Hence this Pointer<DistributionImplementation> 859 mutable Implementation p_standardDistribution_; 860 861 /** Optimization for the generating function evaluation */ 862 mutable Bool isAlreadyCreatedGeneratingFunction_; 863 mutable UniVariatePolynomial generatingFunction_; 864 865 #ifndef SWIG 866 867 // Class used to wrap the computePDF() method for interpolation purpose 868 class PDFWrapper: public EvaluationImplementation 869 { 870 public: PDFWrapper(const DistributionImplementation * p_distribution)871 PDFWrapper(const DistributionImplementation * p_distribution) 872 : EvaluationImplementation() 873 , p_distribution_(p_distribution) 874 { 875 // Nothing to do 876 } 877 clone() const878 PDFWrapper * clone() const override 879 { 880 return new PDFWrapper(*this); 881 } 882 operator ()(const Point & point) const883 Point operator() (const Point & point) const override 884 { 885 return Point(1, p_distribution_->computePDF(point)); 886 } 887 operator ()(const Sample & sample) const888 Sample operator() (const Sample & sample) const override 889 { 890 return p_distribution_->computePDF(sample); 891 }; 892 getInputDimension() const893 UnsignedInteger getInputDimension() const override 894 { 895 return p_distribution_->getDimension(); 896 } 897 getOutputDimension() const898 UnsignedInteger getOutputDimension() const override 899 { 900 return 1; 901 } 902 getInputDescription() const903 Description getInputDescription() const override 904 { 905 return p_distribution_->getDescription(); 906 } 907 getOutputDescription() const908 Description getOutputDescription() const override 909 { 910 return Description(1, "pdf"); 911 } 912 __repr__() const913 String __repr__() const override 914 { 915 OSS oss; 916 oss << "PDFWrapper(" << p_distribution_->__str__() << ")"; 917 return oss; 918 } 919 __str__(const String &) const920 String __str__(const String & ) const override 921 { 922 OSS oss; 923 oss << "PDFWrapper(" << p_distribution_->__str__() << ")"; 924 return oss; 925 } 926 927 private: 928 const DistributionImplementation * p_distribution_; 929 }; // class PDFWrapper 930 931 // Class used to wrap the computeLogPDF() method for interpolation purpose 932 class LogPDFWrapper: public EvaluationImplementation 933 { 934 public: LogPDFWrapper(const DistributionImplementation * p_distribution)935 LogPDFWrapper(const DistributionImplementation * p_distribution) 936 : EvaluationImplementation() 937 , p_distribution_(p_distribution) 938 { 939 // Nothing to do 940 } 941 clone() const942 LogPDFWrapper * clone() const override 943 { 944 return new LogPDFWrapper(*this); 945 } 946 operator ()(const Point & point) const947 Point operator() (const Point & point) const override 948 { 949 return Point(1, p_distribution_->computeLogPDF(point)); 950 } 951 operator ()(const Sample & sample) const952 Sample operator() (const Sample & sample) const override 953 { 954 return p_distribution_->computeLogPDF(sample); 955 }; 956 getInputDimension() const957 UnsignedInteger getInputDimension() const override 958 { 959 return p_distribution_->getDimension(); 960 } 961 getOutputDimension() const962 UnsignedInteger getOutputDimension() const override 963 { 964 return 1; 965 } 966 getInputDescription() const967 Description getInputDescription() const override 968 { 969 return p_distribution_->getDescription(); 970 } 971 getOutputDescription() const972 Description getOutputDescription() const override 973 { 974 return Description(1, "logpdf"); 975 } 976 __repr__() const977 String __repr__() const override 978 { 979 OSS oss; 980 oss << "LogPDFWrapper(" << p_distribution_->__str__() << ")"; 981 return oss; 982 } 983 __str__(const String &) const984 String __str__(const String & ) const override 985 { 986 OSS oss; 987 oss << "LogPDFWrapper(" << p_distribution_->__str__() << ")"; 988 return oss; 989 } 990 991 private: 992 const DistributionImplementation * p_distribution_; 993 }; // class LogPDFWrapper 994 995 // Structure used to wrap the computeCDF() method for interpolation purpose 996 class CDFWrapper: public EvaluationImplementation 997 { 998 public: CDFWrapper(const DistributionImplementation * p_distribution)999 CDFWrapper(const DistributionImplementation * p_distribution) 1000 : EvaluationImplementation() 1001 , p_distribution_(p_distribution) 1002 { 1003 // Nothing to do 1004 } 1005 clone() const1006 CDFWrapper * clone() const override 1007 { 1008 return new CDFWrapper(*this); 1009 } 1010 operator ()(const Point & point) const1011 Point operator() (const Point & point) const override 1012 { 1013 return computeCDF(point); 1014 } 1015 operator ()(const Sample & sample) const1016 Sample operator() (const Sample & sample) const override 1017 { 1018 return computeCDF(sample); 1019 } 1020 computeCDF(const Point & point) const1021 Point computeCDF(const Point & point) const 1022 { 1023 return Point(1, p_distribution_->computeCDF(point)); 1024 }; 1025 computeCDF(const Sample & sample) const1026 Sample computeCDF(const Sample & sample) const 1027 { 1028 return p_distribution_->computeCDF(sample); 1029 }; 1030 getInputDimension() const1031 UnsignedInteger getInputDimension() const override 1032 { 1033 return p_distribution_->getDimension(); 1034 } 1035 getOutputDimension() const1036 UnsignedInteger getOutputDimension() const override 1037 { 1038 return 1; 1039 } 1040 getInputDescription() const1041 Description getInputDescription() const override 1042 { 1043 return p_distribution_->getDescription(); 1044 } 1045 getOutputDescription() const1046 Description getOutputDescription() const override 1047 { 1048 return Description(1, "cdf"); 1049 } 1050 __repr__() const1051 String __repr__() const override 1052 { 1053 OSS oss; 1054 oss << "CDFWrapper(" << p_distribution_->__str__() << ")"; 1055 return oss; 1056 } 1057 __str__(const String &) const1058 String __str__(const String & ) const override 1059 { 1060 OSS oss; 1061 oss << "CDFWrapper(" << p_distribution_->__str__() << ")"; 1062 return oss; 1063 } 1064 1065 private: 1066 const DistributionImplementation * p_distribution_; 1067 }; // class CDFWrapper 1068 1069 // Structure used to implement the computeQuantile() method efficiently 1070 class QuantileWrapper: public EvaluationImplementation 1071 { 1072 public: QuantileWrapper(const Collection<Pointer<DistributionImplementation>> marginals,const DistributionImplementation * p_distribution)1073 QuantileWrapper(const Collection< Pointer<DistributionImplementation> > marginals, 1074 const DistributionImplementation * p_distribution) 1075 : EvaluationImplementation() 1076 , marginals_(marginals) 1077 , p_distribution_(p_distribution) 1078 , dimension_(p_distribution->getDimension()) 1079 { 1080 // Nothing to do 1081 } 1082 clone() const1083 QuantileWrapper * clone() const override 1084 { 1085 return new QuantileWrapper(*this); 1086 } 1087 operator ()(const Point & point) const1088 Point operator() (const Point & point) const override 1089 { 1090 return p_distribution_->computeQuantile(point[0]); 1091 } 1092 operator ()(const Sample & sample) const1093 Sample operator() (const Sample & sample) const override 1094 { 1095 return p_distribution_->computeQuantile(sample.asPoint()); 1096 } 1097 computeDiagonal(const Point & u) const1098 Point computeDiagonal(const Point & u) const 1099 { 1100 const Scalar cdf = p_distribution_->computeCDF(diagonalToSpace(u[0])); 1101 LOGDEBUG(OSS(false) << "in DistributionImplementation::QuantileWrapper::computeDiagonal, u=" << u << ", cdf=" << cdf); 1102 return Point(1, cdf); 1103 } 1104 diagonalToSpace(const Scalar tau) const1105 Point diagonalToSpace(const Scalar tau) const 1106 { 1107 Point x(dimension_); 1108 for (UnsignedInteger i = 0; i < dimension_; ++i) x[i] = marginals_[i]->computeQuantile(tau)[0]; 1109 LOGDEBUG(OSS(false) << "in DistributionImplementation::QuantileWrapper::diagonalToSpace, tau=" << tau << ", x=" << x); 1110 return x; 1111 } 1112 getInputDimension() const1113 UnsignedInteger getInputDimension() const override 1114 { 1115 return 1; 1116 } 1117 getOutputDimension() const1118 UnsignedInteger getOutputDimension() const override 1119 { 1120 return p_distribution_->getDimension(); 1121 } 1122 getInputDescription() const1123 Description getInputDescription() const override 1124 { 1125 return Description(1, "quantile"); 1126 } 1127 getOutputDescription() const1128 Description getOutputDescription() const override 1129 { 1130 return p_distribution_->getDescription(); 1131 } 1132 __repr__() const1133 String __repr__() const override 1134 { 1135 OSS oss; 1136 oss << "QuantileWrapper(" << p_distribution_->__str__() << ")"; 1137 return oss; 1138 } 1139 __str__(const String &) const1140 String __str__(const String & ) const override 1141 { 1142 OSS oss; 1143 oss << "QuantileWrapper(" << p_distribution_->__str__() << ")"; 1144 return oss; 1145 } 1146 1147 private: 1148 const Collection< Pointer<DistributionImplementation> > marginals_; 1149 const DistributionImplementation * p_distribution_; 1150 const UnsignedInteger dimension_; 1151 }; // struct QuantileWrapper 1152 1153 // Structure used to implement the computeInverseSurvivalFunction() method efficiently 1154 class SurvivalFunctionWrapper: public EvaluationImplementation 1155 { 1156 public: SurvivalFunctionWrapper(const Collection<Pointer<DistributionImplementation>> marginals,const DistributionImplementation * p_distribution)1157 SurvivalFunctionWrapper(const Collection< Pointer<DistributionImplementation> > marginals, 1158 const DistributionImplementation * p_distribution) 1159 : EvaluationImplementation() 1160 , marginals_(marginals) 1161 , p_distribution_(p_distribution) 1162 , dimension_(p_distribution->getDimension()) 1163 { 1164 // Nothing to do 1165 } 1166 clone() const1167 SurvivalFunctionWrapper * clone() const override 1168 { 1169 return new SurvivalFunctionWrapper(*this); 1170 } 1171 operator ()(const Point & point) const1172 Point operator() (const Point & point) const override 1173 { 1174 return Point(1, p_distribution_->computeSurvivalFunction(point)); 1175 } 1176 operator ()(const Sample & sample) const1177 Sample operator() (const Sample & sample) const override 1178 { 1179 return p_distribution_->computeSurvivalFunction(sample); 1180 } 1181 computeDiagonal(const Point & u) const1182 Point computeDiagonal(const Point & u) const 1183 { 1184 const Scalar survival = p_distribution_->computeSurvivalFunction(diagonalToSpace(u[0])); 1185 LOGDEBUG(OSS(false) << "in DistributionImplementation::InverseSurvivalFunctionWrapper::computeDiagonal, u=" << u << ", survival=" << survival); 1186 return Point(1, survival); 1187 } 1188 diagonalToSpace(const Scalar tau) const1189 Point diagonalToSpace(const Scalar tau) const 1190 { 1191 Point x(dimension_); 1192 for (UnsignedInteger i = 0; i < dimension_; ++i) x[i] = marginals_[i]->computeQuantile(tau, true)[0]; 1193 LOGDEBUG(OSS(false) << "in DistributionImplementation::InverseSurvivalFunctionWrapper::diagonalToSpace, tau=" << tau << ", x=" << x); 1194 return x; 1195 } 1196 getInputDimension() const1197 UnsignedInteger getInputDimension() const override 1198 { 1199 return p_distribution_->getDimension(); 1200 } 1201 getOutputDimension() const1202 UnsignedInteger getOutputDimension() const override 1203 { 1204 return 1; 1205 } 1206 getInputDescription() const1207 Description getInputDescription() const override 1208 { 1209 return p_distribution_->getDescription(); 1210 } 1211 getOutputDescription() const1212 Description getOutputDescription() const override 1213 { 1214 return Description(1, "survival function"); 1215 } 1216 __repr__() const1217 String __repr__() const override 1218 { 1219 OSS oss; 1220 oss << "SurvivalFunctionWrapper(" << p_distribution_->__str__() << ")"; 1221 return oss; 1222 } 1223 __str__(const String &) const1224 String __str__(const String & ) const override 1225 { 1226 OSS oss; 1227 oss << "SurvivalFunctionWrapper(" << p_distribution_->__str__() << ")"; 1228 return oss; 1229 } 1230 1231 private: 1232 const Collection< Pointer<DistributionImplementation> > marginals_; 1233 const DistributionImplementation * p_distribution_; 1234 const UnsignedInteger dimension_; 1235 }; // struct SurvivalFunctionWrapper 1236 1237 class MinimumVolumeLevelSetEvaluation: public EvaluationImplementation 1238 { 1239 public: 1240 // Here we use a smart pointer instead of a const C++ pointer because the life-cycle of the 1241 // object goes outside of the calling method MinimumVolumeLevelSetEvaluation(const DistributionImplementation::Implementation & p_distribution)1242 MinimumVolumeLevelSetEvaluation(const DistributionImplementation::Implementation & p_distribution) 1243 : EvaluationImplementation() 1244 , p_distribution_(p_distribution) 1245 { 1246 // Nothing to do 1247 } 1248 clone() const1249 MinimumVolumeLevelSetEvaluation * clone() const override 1250 { 1251 return new MinimumVolumeLevelSetEvaluation(*this); 1252 } 1253 1254 // The minimum volume level A(p) set is such that A(p)={x\in R^n | y(x) <= y_p} 1255 // where y(x)=-\log X and y_p is the p-quantile of Y=pdf(X) operator ()(const Point & point) const1256 Point operator() (const Point & point) const override 1257 { 1258 const Scalar value = -p_distribution_->computeLogPDF(point); 1259 return Point(1, value); 1260 } 1261 operator ()(const Sample & sample) const1262 Sample operator() (const Sample & sample) const override 1263 { 1264 return p_distribution_->computeLogPDF(sample) * (-1.0); 1265 } 1266 getInputDimension() const1267 UnsignedInteger getInputDimension() const override 1268 { 1269 return p_distribution_->getDimension(); 1270 } 1271 getOutputDimension() const1272 UnsignedInteger getOutputDimension() const override 1273 { 1274 return 1; 1275 } 1276 getInputDescription() const1277 Description getInputDescription() const override 1278 { 1279 return p_distribution_->getDescription(); 1280 } 1281 getOutputDescription() const1282 Description getOutputDescription() const override 1283 { 1284 return Description(1, "-logPDF"); 1285 } 1286 getDescription() const1287 Description getDescription() const override 1288 { 1289 Description description(getInputDescription()); 1290 description.add(getOutputDescription()); 1291 return description; 1292 } 1293 __repr__() const1294 String __repr__() const override 1295 { 1296 OSS oss; 1297 oss << "MinimumVolumeLevelSetEvaluation(" << p_distribution_->__str__() << ")"; 1298 return oss; 1299 } 1300 __str__(const String &) const1301 String __str__(const String & ) const override 1302 { 1303 OSS oss; 1304 oss << "MinimumVolumeLevelSetEvaluation(" << p_distribution_->__str__() << ")"; 1305 return oss; 1306 } 1307 1308 private: 1309 const DistributionImplementation::Implementation p_distribution_; 1310 }; // class MinimumVolumeLevelSetEvaluation 1311 1312 class MinimumVolumeLevelSetGradient: public GradientImplementation 1313 { 1314 public: 1315 // Here we use a smart pointer instead of a const C++ pointer because the life-cycle of the 1316 // object goes outside of the calling method MinimumVolumeLevelSetGradient(const DistributionImplementation::Implementation & p_distribution)1317 MinimumVolumeLevelSetGradient(const DistributionImplementation::Implementation & p_distribution) 1318 : GradientImplementation() 1319 , p_distribution_(p_distribution) 1320 { 1321 // Nothing to do 1322 } 1323 clone() const1324 MinimumVolumeLevelSetGradient * clone() const override 1325 { 1326 return new MinimumVolumeLevelSetGradient(*this); 1327 } 1328 gradient(const Point & point) const1329 Matrix gradient(const Point & point) const override 1330 { 1331 const Scalar pdf = p_distribution_->computePDF(point); 1332 if (pdf == 0) return Matrix(getInputDimension(), getOutputDimension()); 1333 const Point value = p_distribution_->computeDDF(point) * (-1.0 / pdf); 1334 return MatrixImplementation(getInputDimension(), getOutputDimension(), value); 1335 } 1336 getInputDimension() const1337 UnsignedInteger getInputDimension() const override 1338 { 1339 return p_distribution_->getDimension(); 1340 } 1341 getOutputDimension() const1342 UnsignedInteger getOutputDimension() const override 1343 { 1344 return 1; 1345 } 1346 getInputDescription() const1347 Description getInputDescription() const 1348 { 1349 return p_distribution_->getDescription(); 1350 } 1351 getOutputDescription() const1352 Description getOutputDescription() const 1353 { 1354 return Description(1, "-logPDF"); 1355 } 1356 getDescription() const1357 Description getDescription() const 1358 { 1359 Description description(getInputDescription()); 1360 description.add(getOutputDescription()); 1361 return description; 1362 } 1363 __repr__() const1364 String __repr__() const override 1365 { 1366 OSS oss; 1367 oss << "MinimumVolumeLevelSetGradient(" << p_distribution_->__str__() << ")"; 1368 return oss; 1369 } 1370 __str__(const String &) const1371 String __str__(const String & ) const override 1372 { 1373 OSS oss; 1374 oss << "MinimumVolumeLevelSetGradient(" << p_distribution_->__str__() << ")"; 1375 return oss; 1376 } 1377 1378 private: 1379 const DistributionImplementation::Implementation p_distribution_; 1380 }; // class MinimumVolumeLevelSetGradient 1381 1382 class CovarianceWrapper: public EvaluationImplementation 1383 { 1384 public: CovarianceWrapper(const DistributionImplementation::Implementation & p_distribution,const Scalar muI,const Scalar muJ)1385 CovarianceWrapper(const DistributionImplementation::Implementation & p_distribution, 1386 const Scalar muI, 1387 const Scalar muJ) 1388 : EvaluationImplementation() 1389 , p_distribution_(p_distribution) 1390 , muI_(muI) 1391 , muJ_(muJ) 1392 { 1393 // Nothing to do 1394 } 1395 clone() const1396 CovarianceWrapper * clone() const override 1397 { 1398 return new CovarianceWrapper(*this); 1399 } 1400 operator ()(const Point & point) const1401 Point operator() (const Point & point) const override 1402 { 1403 return Point(1, (point[0] - muI_) * (point[1] - muJ_) * p_distribution_->computePDF(point)); 1404 } 1405 operator ()(const Sample & sample) const1406 Sample operator() (const Sample & sample) const override 1407 { 1408 const UnsignedInteger size = sample.getSize(); 1409 Sample result(p_distribution_->computePDF(sample)); 1410 for (UnsignedInteger i = 0; i < size; ++i) result(i, 0) *= (sample(i, 0) - muI_) * (sample(i, 1) - muJ_); 1411 return result; 1412 } 1413 getInputDimension() const1414 UnsignedInteger getInputDimension() const override 1415 { 1416 return 2; 1417 } 1418 getOutputDimension() const1419 UnsignedInteger getOutputDimension() const override 1420 { 1421 return 1; 1422 } 1423 getInputDescription() const1424 Description getInputDescription() const override 1425 { 1426 return Description::BuildDefault(2, "x"); 1427 } 1428 getOutputDescription() const1429 Description getOutputDescription() const override 1430 { 1431 return Description(1, "c"); 1432 } 1433 __repr__() const1434 String __repr__() const override 1435 { 1436 OSS oss; 1437 oss << "CovarianceWrapper(" << p_distribution_->__str__() << ")"; 1438 return oss; 1439 } 1440 __str__(const String &) const1441 String __str__(const String & ) const override 1442 { 1443 OSS oss; 1444 oss << "CovarianceWrapper(" << p_distribution_->__str__() << ")"; 1445 return oss; 1446 } 1447 1448 private: 1449 const DistributionImplementation::Implementation p_distribution_; 1450 const Scalar muI_; 1451 const Scalar muJ_; 1452 }; // class CovarianceWrapper 1453 1454 // Class used to wrap the computeConditionalPDF() method for the computation of the conditional CDF 1455 class ShiftedMomentWrapper: public EvaluationImplementation 1456 { 1457 public: ShiftedMomentWrapper(const UnsignedInteger n,const Scalar shift,const DistributionImplementation::Implementation & p_distribution)1458 ShiftedMomentWrapper(const UnsignedInteger n, 1459 const Scalar shift, 1460 const DistributionImplementation::Implementation & p_distribution) 1461 : EvaluationImplementation() 1462 , n_(1.0 * n) 1463 , shift_(shift) 1464 , p_distribution_(p_distribution) 1465 { 1466 // Nothing to do 1467 }; 1468 clone() const1469 ShiftedMomentWrapper * clone() const override 1470 { 1471 return new ShiftedMomentWrapper(*this); 1472 } 1473 operator ()(const Point & point) const1474 Point operator() (const Point & point) const override 1475 { 1476 const Scalar power = std::pow(point[0] - shift_, n_); 1477 const Scalar pdf = p_distribution_->computePDF(point); 1478 return Point(1, power * pdf); 1479 }; 1480 operator ()(const Sample & sample) const1481 Sample operator() (const Sample & sample) const override 1482 { 1483 const UnsignedInteger size = sample.getSize(); 1484 Sample result(size, 1); 1485 const Sample pdf(p_distribution_->computePDF(sample)); 1486 for (UnsignedInteger i = 0; i < size; ++i) 1487 result(i, 0) = std::pow(sample(i, 0) - shift_, n_) * pdf(i, 0); 1488 return result; 1489 }; 1490 getInputDimension() const1491 UnsignedInteger getInputDimension() const override 1492 { 1493 return 1; 1494 } 1495 getOutputDimension() const1496 UnsignedInteger getOutputDimension() const override 1497 { 1498 return 1; 1499 } 1500 __repr__() const1501 String __repr__() const override 1502 { 1503 OSS oss; 1504 oss << "ShiftedMomentWrapper(" << p_distribution_->__str__() << ")"; 1505 return oss; 1506 } 1507 __str__(const String &) const1508 String __str__(const String & ) const override 1509 { 1510 OSS oss; 1511 oss << "ShiftedMomentWrapper(" << p_distribution_->__str__() << ")"; 1512 return oss; 1513 } 1514 1515 private: 1516 const Scalar n_; 1517 const Scalar shift_; 1518 const DistributionImplementation::Implementation & p_distribution_; 1519 }; // class ShiftedMomentWrapper 1520 1521 // Class used to wrap the computeConditionalPDF() method for the computation of the conditional CDF 1522 class ConditionalPDFWrapper: public UniVariateFunctionImplementation 1523 { 1524 public: ConditionalPDFWrapper(const DistributionImplementation::Implementation p_distribution)1525 ConditionalPDFWrapper(const DistributionImplementation::Implementation p_distribution) 1526 : UniVariateFunctionImplementation() 1527 , y_(1, 0.0) 1528 , p_distribution_(p_distribution) 1529 { 1530 // Nothing to do 1531 }; 1532 operator ()(const Scalar x) const1533 Scalar operator() (const Scalar x) const override 1534 { 1535 Collection<Scalar> z(y_); 1536 z.add(x); 1537 return p_distribution_->computePDF(z); 1538 }; 1539 setParameter(const Point & parameters)1540 void setParameter(const Point & parameters) 1541 { 1542 y_ = parameters; 1543 } 1544 getParameter() const1545 Point getParameter() const 1546 { 1547 return y_; 1548 } 1549 __repr__() const1550 String __repr__() const override 1551 { 1552 OSS oss; 1553 oss << "ConditionalPDFWrapper(" << p_distribution_->__str__() << ")"; 1554 return oss; 1555 } 1556 __str__(const String &) const1557 String __str__(const String & ) const override 1558 { 1559 OSS oss; 1560 oss << "ConditionalPDFWrapper(" << p_distribution_->__str__() << ")"; 1561 return oss; 1562 } 1563 1564 private: 1565 Collection<Scalar> y_; 1566 const DistributionImplementation::Implementation p_distribution_; 1567 }; // class ConditionalPDFWrapper 1568 1569 // Class used to wrap the computeConditionalCDF() method for the computation of the conditional quantile 1570 class ConditionalCDFWrapper: public UniVariateFunctionImplementation 1571 { 1572 public: ConditionalCDFWrapper(const DistributionImplementation * p_distribution)1573 ConditionalCDFWrapper(const DistributionImplementation * p_distribution) 1574 : UniVariateFunctionImplementation() 1575 , y_(1, 0.0) 1576 , p_distribution_(p_distribution) 1577 { 1578 // Nothing to do 1579 }; 1580 clone() const1581 ConditionalCDFWrapper * clone() const override 1582 { 1583 return new ConditionalCDFWrapper(*this); 1584 } 1585 operator ()(const Scalar x) const1586 Scalar operator() (const Scalar x) const override 1587 { 1588 return p_distribution_->computeConditionalCDF(x, y_); 1589 }; 1590 setParameter(const Point & parameters)1591 void setParameter(const Point & parameters) 1592 { 1593 y_ = parameters; 1594 } 1595 getParameter() const1596 Point getParameter() const 1597 { 1598 return y_; 1599 } 1600 __repr__() const1601 String __repr__() const override 1602 { 1603 OSS oss; 1604 oss << "ConditionalCDFWrapper(" << p_distribution_->__str__() << ")"; 1605 return oss; 1606 } 1607 __str__(const String &) const1608 String __str__(const String & ) const override 1609 { 1610 OSS oss; 1611 oss << "ConditionalCDFWrapper(" << p_distribution_->__str__() << ")"; 1612 return oss; 1613 } 1614 1615 private: 1616 Point y_; 1617 const DistributionImplementation * p_distribution_; 1618 }; // class ConditionalCDFWrapper 1619 1620 // Class used to compute entropy 1621 class EntropyKernel: public EvaluationImplementation 1622 { 1623 public: EntropyKernel(const DistributionImplementation * p_distribution)1624 EntropyKernel(const DistributionImplementation * p_distribution) 1625 : EvaluationImplementation() 1626 , p_distribution_(p_distribution) 1627 { 1628 // Nothing to do 1629 } 1630 clone() const1631 EntropyKernel * clone() const override 1632 { 1633 return new EntropyKernel(*this); 1634 } 1635 operator ()(const Point & point) const1636 Point operator() (const Point & point) const override 1637 { 1638 const Scalar logPDF = p_distribution_->computeLogPDF(point); 1639 return Point(1, -std::exp(logPDF) * logPDF); 1640 } 1641 operator ()(const Sample & sample) const1642 Sample operator() (const Sample & sample) const override 1643 { 1644 const UnsignedInteger size = sample.getSize(); 1645 const Point logPDF(p_distribution_->computeLogPDF(sample).asPoint()); 1646 Sample result(size, 1); 1647 for (UnsignedInteger i = 0; i < size; ++i) 1648 { 1649 const Scalar logPDFI = logPDF[i]; 1650 result(i, 0) = -std::exp(logPDFI) * logPDFI; 1651 } 1652 return result; 1653 }; 1654 getInputDimension() const1655 UnsignedInteger getInputDimension() const override 1656 { 1657 return p_distribution_->getDimension(); 1658 } 1659 getOutputDimension() const1660 UnsignedInteger getOutputDimension() const override 1661 { 1662 return 1; 1663 } 1664 getInputDescription() const1665 Description getInputDescription() const override 1666 { 1667 return p_distribution_->getDescription(); 1668 } 1669 getOutputDescription() const1670 Description getOutputDescription() const override 1671 { 1672 return Description(1, "entropyKernel"); 1673 } 1674 __repr__() const1675 String __repr__() const override 1676 { 1677 OSS oss; 1678 oss << "EntropyKernel(" << p_distribution_->__str__() << ")"; 1679 return oss; 1680 } 1681 __str__(const String &) const1682 String __str__(const String & ) const override 1683 { 1684 OSS oss; 1685 oss << "EntropyKernel(" << p_distribution_->__str__() << ")"; 1686 return oss; 1687 } 1688 1689 private: 1690 const DistributionImplementation * p_distribution_; 1691 }; // class EntropyKernel 1692 1693 // Class used to wrap the computePDF()**2 method for integration purpose 1694 class PDFSquaredWrapper : public EvaluationImplementation 1695 { 1696 public: PDFSquaredWrapper(const DistributionImplementation * p_distribution)1697 PDFSquaredWrapper(const DistributionImplementation *p_distribution) 1698 : EvaluationImplementation(), p_distribution_(p_distribution) 1699 { 1700 // Nothing to do 1701 } 1702 clone() const1703 PDFSquaredWrapper *clone() const override 1704 { 1705 return new PDFSquaredWrapper(*this); 1706 } 1707 operator ()(const Point & point) const1708 Point operator()(const Point &point) const override 1709 { 1710 const Scalar pdf = p_distribution_->computePDF(point); 1711 return Point(1, pdf * pdf); 1712 } 1713 operator ()(const Sample & sample) const1714 Sample operator()(const Sample &sample) const override 1715 { 1716 // Keep the parallelism of computePDF 1717 Sample pdfSquared(p_distribution_->computePDF(sample)); 1718 for (UnsignedInteger i = 0; i < pdfSquared.getSize(); ++i) 1719 pdfSquared(i, 0) = pdfSquared(i, 0) * pdfSquared(i, 0); 1720 return pdfSquared; 1721 } 1722 getInputDimension() const1723 UnsignedInteger getInputDimension() const override 1724 { 1725 return p_distribution_->getDimension(); 1726 } 1727 getOutputDimension() const1728 UnsignedInteger getOutputDimension() const override 1729 { 1730 return 1; 1731 } 1732 getInputDescription() const1733 Description getInputDescription() const override 1734 { 1735 return p_distribution_->getDescription(); 1736 } 1737 getOutputDescription() const1738 Description getOutputDescription() const override 1739 { 1740 return Description(1, "pdfSquared"); 1741 } 1742 __repr__() const1743 String __repr__() const override 1744 { 1745 OSS oss; 1746 oss << "PDFSquaredWrapper(" << p_distribution_->__str__() << ")"; 1747 return oss; 1748 } 1749 __str__(const String &) const1750 String __str__(const String &) const override 1751 { 1752 OSS oss; 1753 oss << "PDFSquaredWrapper(" << p_distribution_->__str__() << ")"; 1754 return oss; 1755 } 1756 1757 private: 1758 const DistributionImplementation *p_distribution_; 1759 }; // class PDFSquaredWrapper 1760 1761 #endif 1762 1763 /** The dimension of the distribution */ 1764 UnsignedInteger dimension_; 1765 1766 /** The weight used ONLY by Mixture */ 1767 Scalar weight_; 1768 1769 /** Range of the distribution */ 1770 Interval range_; 1771 1772 /** Description of each component */ 1773 Description description_; 1774 1775 /** Use parallelism */ 1776 Bool isParallel_; 1777 1778 /** Flag to tell if the distribution is indeed a copula */ 1779 Bool isCopula_; 1780 1781 }; /* class DistributionImplementation */ 1782 1783 #ifndef SWIG 1784 1785 /** Product operator */ 1786 OT_API Distribution operator * (const Scalar, 1787 const DistributionImplementation & distribution); 1788 1789 OT_API Distribution operator * (const Scalar, 1790 const DistributionImplementation::Implementation & p_distribution); 1791 1792 /** Division operator */ 1793 OT_API Distribution operator / (const Scalar, 1794 const DistributionImplementation & distribution); 1795 1796 OT_API Distribution operator / (const Scalar, 1797 const DistributionImplementation::Implementation & p_distribution); 1798 1799 /** Addition operator */ 1800 OT_API Distribution operator + (const Scalar, 1801 const DistributionImplementation & distribution); 1802 1803 OT_API Distribution operator + (const Scalar, 1804 const DistributionImplementation::Implementation & p_distribution); 1805 1806 /** Subtraction operator */ 1807 OT_API Distribution operator - (const Scalar, 1808 const DistributionImplementation & distribution); 1809 1810 OT_API Distribution operator - (const Scalar, 1811 const DistributionImplementation::Implementation & p_distribution); 1812 1813 OT_API Distribution operator - (const DistributionImplementation & distribution); 1814 1815 OT_API Distribution operator - (const DistributionImplementation::Implementation & p_distribution); 1816 1817 1818 #endif 1819 1820 OT_API Distribution maximum(const DistributionImplementation::Implementation & p_left, 1821 const DistributionImplementation::Implementation & p_right); 1822 1823 OT_API Distribution maximum(const DistributionImplementation & left, 1824 const DistributionImplementation::Implementation & p_right); 1825 1826 OT_API Distribution maximum(const DistributionImplementation::Implementation & p_left, 1827 const DistributionImplementation & right); 1828 1829 OT_API Distribution maximum(const DistributionImplementation & left, 1830 const DistributionImplementation & right); 1831 1832 1833 END_NAMESPACE_OPENTURNS 1834 1835 #endif /* OPENTURNS_DISTRIBUTIONIMPLEMENTATION_HXX */ 1836