R/Kernel_Tricube.R

#' @title Tricube Kernel
#'
#' @description Mathematical and statistical functions for the Tricube kernel defined by the pdf,
#' \deqn{f(x) = 70/81(1 - |x|^3)^3}
#' over the support \eqn{x \in (-1,1)}{x \epsilon (-1,1)}.
#'
#' @details The quantile function is omitted as no closed form analytic expressions could
#' be found, decorate with FunctionImputation for numeric results.
#'
#' @name Tricube
#' @template class_distribution
#' @template class_kernel
#' @template method_pdfsquared2Norm
#'
#' @export
Tricube <- R6Class("Tricube",
  inherit = Kernel, lock_objects = F,
  public = list(
    name = "Tricube",
    short_name = "Tric",
    description = "Tricube Kernel",
    #' @description
    #' The squared 2-norm of the pdf is defined by
    #' \deqn{\int_a^b (f_X(u))^2 du}
    #' where X is the Distribution, \eqn{f_X} is its pdf and \eqn{a, b}
    #' are the distribution support limits.
    pdfSquared2Norm = function(x = 0, upper = Inf) {

      fixed1 <- function(x, upper) {

        ret1 <- (upper * (- (4900 * x^9) / 6561 + (4900 * x^6) / 2187 -
                            (4900 * x^3) / 2187 + 4900 / 6561) +
                   upper^2 * ((2450 * x^8) / 729 - (4900 * x^5) / 729 + (2450 * x^2) / 729) +
                   upper^3 * (- (19600 * x^7) / 2187 + (24500 * x ^4) / 2187 - (4900 * x) / 2187) +
                   upper^4 * (- (1225 * x^9) / 2187 + (37975 * x^6) / 2187 -
                               (28175 * x^3) / 2187 + 2450 / 2187) +
                   upper^5 * ((980 * x^8) / 243 - (19600 * x^5) / 729 + (7840 * x^2) / 729) +
                   upper^6 * (- (9800 * x^7) / 729 + (71050 * x^4) / 2187 - (12250 * x) / 2187) +
                   upper^7 * (- (700 * x^9) / 2187 + (20300 * x^6) / 729 -
                                (63700 * x^3) / 2187 + 3500 / 2187) +
                   upper^8 * ((1225 * x^8) / 486 - (9800 * x^5) / 243 + (13475 * x^2) / 729) +
                   upper^9 * (- (19600 * x^7) / 2187 + (93100 * x^4) / 2187 - (49000 * x) / 6561) +
                   upper^10 * (- (490 * x^9) / 6561 + (41650 * x^6) / 2187 -
                                 (71050 * x^3) / 2187 + 9800 / 6561) +
                   upper^11 * ((4900 * x^8) / 8019 - (19600 * x^5) / 729 + (137200 * x^2) / 8019) +
                   upper^12 * (- (4900 * x^7) / 2187 + (57575 * x^4) / 2187 - (12250 * x) / 2187) +
                   upper^13 * ((137200 * x^6) / 28431 - (39200 * x^3) / 2187 + 24500 / 28431) +
                   upper^14 * ((5950 * x^2) / 729 - (4900 * x^5) / 729) +
                   upper^15 * ((13720 * x^4) / 2187 - (4900 * x) / 2187) +
                   upper^16 * (1225 / 4374 - (8575 * x^3) / 2187) +
                   (19600 * upper^17 * x^2) / 12393 - (2450 * upper^18 * x) / 6561 +
                   (4900 * upper^19) / 124659)
        return(ret1)

      }

      fixed2 <- function(x, upper) {

        ret2 <- abs(upper)^4 * ((1225 * abs(x)^9) / 2187 +
          (30625 * abs(x)^6) / 2187 - (20825 * abs(x)^3) / 2187) +
          abs(upper)^10 * ((490 * abs(x)^9) / 6561 +
          (40670 * abs(x)^6) / 2187 + (12250 * abs(x)^3) / 2187) +
          abs(upper)^7 * (- (700 * abs(x)^9) / 2187 -
          (700 * abs(x)^6) / 27 + (20300 * abs(x)^3) / 2187 - 700 / 2187) +
          abs(upper) * (- (4900 * abs(x)^9) / 6561 +
          (4900 * abs(x)^6) / 2187 - (4900 * abs(x)^3) / 2187 + 4900 / 6561) +
          abs(upper)^2 * ((2450 * abs(x)^8) / 729 -
          (4900 * abs(x)^5) / 729 + (2450 * abs(x)^2) / 729) +
          abs(upper)^8 * ((1225 * abs(x)^8) / 486 + (2450 * abs(x)^5) / 81 -
          (4900 * abs(x)^2) / 729) +
          abs(upper)^11 * (- (4900 * abs(x)^8) / 8019 - (196000 * abs(x)^5) / 8019 +
          (9800 * abs(x)^2) / 8019) +
          abs(upper)^5 * (- (980 * abs(x)^8) / 243 - (7840 * abs(x)^5) / 729 +
          (1960 * abs(x)^2) / 729) +
          abs(upper)^6 * ((9800 * abs(x)^7) / 729 - (2450 * abs(x)^4) / 2187 +
          (2450 * abs(x)) / 2187) +
          abs(upper)^12 * ((4900 * abs(x)^7) / 2187 + (45325 * abs(x)^4) / 2187 -
          (2450 * abs(x)) / 2187) +
          abs(upper)^3 * (- (19600 * abs(x)^7) / 2187 + (24500 * abs(x)^4) / 2187 -
                            (4900 * abs(x)) / 2187) +
          abs(upper)^9 * (- (19600 * abs(x)^7) / 2187 - (4900 * abs(x)^4) / 243 +
                            (9800 * abs(x)) / 6561) +
          abs(upper)^13 * (- (137200 * abs(x)^6) / 28431 - (313600 * abs(x)^3) / 28431
                           + 4900 / 28431) +
          abs(upper)^14 * ((4900 * abs(x)^5) / 729 + (2450 * abs(x)^2) / 729) +
          abs(upper)^15 * (- (13720 * abs(x)^4) / 2187 - (980 * abs(x)) / 2187) +
          (8575 * abs(upper)^16 * abs(x)^3) / 2187 - (19600 * abs(upper)^17 * abs(x)^2) / 12393 +
          (2450 * abs(upper)^18 * abs(x)) / 6561 - (4900 * abs(upper)^19) / 124659
        return(ret2)

      }

      fixed3 <- function(x, upper) {

        ret3 <- upper * ((4900 * x^9) / 6561 + (4900 * x^6) / 2187 +
                           (4900 * x^3) / 2187 + 4900 / 6561) +
          upper^7 * ((700 * (x)^9) / 2187 + (20300 * x^6) / 729 +
                       (63700 * x^3) / 2187 + 3500 / 2187) +
          upper^10 * (- (490 * x^9) / 6561 - (41650 * x^6) / 2187 -
                        (71050 * x^3) / 2187 - 9800 / 6561) +
          upper^4 * (- (1225 * x^9) / 2187 - (37975 * x^6) / 2187 -
                       (28175 * x^3) / 2187 - 2450 / 2187) +
          upper^5 * ((980 * x^8) / 243 + (19600 * x^5) / 729 + (7840 * x^2) / 729) +
          upper^11 * ((4900 * x^8) / 8019 + (19600 * x^5) / 729 + (137200 * x^2) / 8019) +
          upper^8 * (- (1225 * x^8) / 486 - (9800 * x^5) / 243 - (13475 * x^2) / 729) +
          upper^2 * (- (2450 * x^8) / 729 - (4900 * x^5) / 729 - (2450 * x^2) / 729) +
          upper^9 * ((19600 * x^7) / 2187 + (93100 * x^4) / 2187 + (49000 * x) / 6561) +
          upper^3 * ((19600 * x^7) / 2187 + (24500 * x^4) / 2187 + (4900 * x) / 2187) +
          upper^12 * (- (4900 * x^7) / 2187 - (57575 * x^4) / 2187 - (12250 * x) / 2187) +
          upper^6 * (- (9800 * x^7) / 729 - (71050 * x^4) / 2187 - (12250 * x) / 2187) +
          upper^13 * ((137200 * x^6) / 28431 + (39200 * x^3) / 2187 + 24500 / 28431) +
          upper^14 * (- (4900 * x^5) / 729 - (5950 * x^2) / 729) +
          upper^15 * ((13720 * x^4) / 2187 + (4900 * x) / 2187) +
          upper^16 * (- (8575 * x^3) / 2187 - 1225 / 4374) +
          (19600 * upper^17 * x^2) / 12393 -
          (2450 * upper^18 * x) / 6561 + (4900 * upper^19) / 124659

        return(ret3)
      }

      fixed4 <- function(x, upper) {

        ret4 <- abs(upper)^4 *
          ((1225 * abs(x)^9) / 2187 + (30625 * abs(x)^6) / 2187 -
             (20825 * abs(x)^3) / 2187) +
          abs(upper)^10 * ((490 * abs(x)^9) / 6561 + (40670 * abs(x)^6) / 2187 +
                             (12250 * abs(x)^3) / 2187) +
          abs(upper)^7 * (- (700 * abs(x)^9) / 2187 - (700 * abs(x)^6) / 27 +
                            (20300 * abs(x)^3) / 2187 - 700 / 2187) +
          abs(upper) * (- (4900 * abs(x)^9) / 6561 + (4900 * abs(x)^6) / 2187 -
                          (4900 * abs(x)^3) / 2187 + 4900 / 6561) +
          abs(upper)^2 * ((2450 * abs(x)^8) / 729 - (4900 * abs(x)^5) / 729 +
                            (2450 * abs(x)^2) / 729) +
          abs(upper)^8 * ((1225 * abs(x)^8) / 486 + (2450 * abs(x)^5) / 81 -
                            (4900 * abs(x)^2) / 729) +
          abs(upper)^11 * (- (4900 * abs(x)^8) / 8019 - (196000 * abs(x)^5) / 8019 +
                             (9800 * abs(x)^2) / 8019) +
          abs(upper)^5 * (- (980 * abs(x)^8) / 243 - (7840 * abs(x)^5) / 729 +
                            (1960 * abs(x)^2) / 729) +
          abs(upper)^6 * ((9800 * abs(x)^7) / 729 - (2450 * abs(x)^4) / 2187 +
                            (2450 * abs(x)) / 2187) +
          abs(upper)^12 * ((4900 * abs(x)^7) / 2187 + (45325 * abs(x)^4) / 2187 -
                             (2450 * abs(x)) / 2187) +
          abs(upper)^3 * (- (19600 * abs(x)^7) / 2187 + (24500 * abs(x)^4) / 2187 -
                            (4900 * abs(x)) / 2187) +
          abs(upper)^9 * (- (19600 * abs(x)^7) / 2187 - (4900 * abs(x)^4) / 243 +
                            (9800 * abs(x)) / 6561) +
          abs(upper)^13 * (- (137200 * abs(x)^6) / 28431 - (313600 * abs(x)^3) / 28431 +
                             4900 / 28431) +
          abs(upper)^14 * ((4900 * abs(x)^5) / 729 + (2450 * abs(x)^2) / 729) +
          abs(upper)^15 * (- (13720 * abs(x)^4) / 2187 - (980 * abs(x)) / 2187) +
          (8575 * abs(upper)^16 * abs(x)^3) / 2187 - (19600 * abs(upper)^17 * abs(x)^2) / 12393 +
          (2450 * abs(upper)^18 * abs(x)) / 6561 - (4900 * abs(upper)^19) / 124659

        return(ret4)
      }

      fixed5 <- function(x) {

        ret5 <- (70 / 81)^2 *
          (6561 / 6916 - (19683 * abs(x)^2) / 13090 + (9 * abs(x)^4) / 5 -
          (729 * abs(x)^6) / 182 +  (747 * abs(x)^7) / 140 -
          (729 * abs(x)^8) / 220 + (81 * abs(x)^9) / 70 - (31 * abs(x)^10) / 140 +
          (111 * abs(x)^13) / 20020 - (3 * abs(x)^16) / 40040 +
          abs(x)^19 / 461890)
        return(ret5)

      }

      fixed6 <- function(x) {

        ret6 <- 22400 / 20007 - (15680 * abs(x)) / 6561 +
          (7840 * x^2) / 1683 - (2800 * abs(x)^3) / 351 +
          (21560 * x^4) / 2187 - (980 * abs(x)^5) / 99 +
          (2870 *  x^6) / 351 - (11305 * abs(x)^7) / 2187 +
          (245 * x^8) / 99 - (70 * abs(x)^9) / 81 +
          (1085 * x^10) / 6561 - (665 * abs(x)^13) / 312741 +
          (35 * x^16) / 625482 - (245 * abs(x)^19) / 303046029

        return(ret6)
      }

      xl <- length(x)
      ul <- length(upper)
      len <- max(xl, ul)

      ret <- numeric(len)
      for (i in seq(len)) {

        xi <- x[ifelse(i %% xl == 0, xl, i %% xl)]
        ui <- upper[ifelse(i %% ul == 0, ul, i %% ul)]

        if (xi >= 0 & xi <= 1) {
          if (ui == Inf | ui > 1) {
            ret[i] <- fixed5(x = xi)
          }  else if (ui > (xi - 1) & ui <= 0) {
            ret[i] <- ((245 * xi^19) / 303046029 + (35 * xi^13) / 34749 +
                         (35 * xi^9) / 81 -
                      (245 * xi^8) / 198 + (3605 * xi^7) / 2187 -
                        (175 * xi^6) / 117 +
                      (2695 * xi^4) / 2187 - (105 * xi^2) / 187 -
                        (2450 * xi) / 6561 + 175 / 494) +
              fixed1(x = xi, upper = ui)
          } else if (ui >= 0 & ui <= xi) {
            ret[i] <- (245 * xi^19) / 303046029 +
              (35 * xi^13) / 34749 + (35 * xi^9) / 81 +
              - (245 * xi^8) / 198 + (3605 * xi^7) / 2187 -
              (175 * xi^6) / 117 + (2695 * xi^4) / 2187 -
              (105 * xi^2) / 187 - (2450 * xi) / 6561 + 175 / 494 +
              fixed2(x = xi, upper = ui)
          } else if (ui >= xi & ui <= 1) {
            ret[i] <- 175 / 494  - (2450 * xi) / 6561 - (105 * xi^2) / 187  +
              (245 * xi^4) / 2187 -
              (175 * xi^6) / 117 + (3815 * xi^7) / 2187 - (245 * xi^8) / 198 +
              (35 * xi^9) / 81 -
              (1085 * xi^10) / 6561 + (1295 * xi^13) / 312741 - (35 * xi^16) / 625482 +
              (245 * xi^19) / 101015343 + fixed3(upper = ui, x = xi)
          } else {
            ret[i] <- 0
            }
        } else if (xi >= 1 & xi < 2) {
          if (ui == Inf | ui > 1) {
            ret[i] <- fixed6(x = xi)
          } else if (ui > (xi - 1) & ui <= 1) {
            ret[i] <- 11200 / 20007 - (3430 * xi) / 2187 + (3920 * xi^2) / 1683  -
              (1400 * xi^3) / 351 + (12005 * xi^4) / 2187 -
              (490 * xi^5) / 99 + (1435 * xi^6) / 351 - (2135 * xi^7) / 729 +
              (245 * xi^8) / 198 - (35 * xi^9) / 81 + (1085 * xi^10) / 6561 -
              (665 * xi^13) / 312741 + (35 * xi^16) / 625482 -
              (245 * xi^19) / 303046029 + fixed4(upper = ui, x = xi)
          } else {
            ret[i] <- 0
          }
        } else if (xi >= - 2 & xi <= - 1) {
          if (ui > - 1 & ui <= xi + 1) {
            ret[i] <- 11200 / 20007 + (5390 * xi) / 6561 + (3920 * xi^2) / 1683  + (
              1400 * xi^3) / 351 + (3185 * xi^4) / 729 + (490 * xi^5) / 99 +
              (1435 * xi^6) / 351 + (4900 * xi^7) / 2187 + (245 * xi^8) / 198 +
              (35 * xi^9) / 81 - fixed4(upper = ui, x = xi)
          } else if (ui == Inf | ui > xi + 1) {
            ret[i] <- fixed6(x = xi)
          } else {
            ret[i] <- 0
            }
        } else if (xi >= -1 & xi <= 0) {
          if (ui > -1 & ui <= xi) {
            ret[i] <- (- (35 * xi^9) / 81 - (245 * xi^8) / 198 -
                         (4900 * xi^7) / 2187 - (175 * xi^6) / 117 +
                       (2695 * xi^4) / 2187 - (105 * xi^2) / 187 -
                         (2450 * xi) / 6561 + 175 / 494) +
              fixed1(x = xi, upper = ui)
          } else if (ui >= xi & ui <= 0) {
            ret[i] <- - (490 * xi^19) / 303046029 - (35 * xi^16) / 625482 -
              (980 * xi^13) / 312741 -
              (1085 * xi^10) / 6561 - (35 * xi^9) / 81 - (245 * xi^8) / 198 -
              (5110 * xi^7) / 2187 - (175 * xi^6) / 117  +
              (245 * xi^4) / 2187 - (105 * xi^2) / 187 - (2450 * xi) / 6561 +
              175 / 494 - fixed2(upper = ui, x = xi)
          } else if (ui >= 0 & ui <= xi + 1) {
            ret[i] <- - (490 * xi^19) / 303046029 - (35 * xi^16) / 625482 -
              (980 * xi^13) / 312741 - (1085 * xi^10) / 6561 -
              (35 * xi^9) / 81 - (245 * xi^8) / 198 -
              (5110 * xi^7) / 2187 - (175 * xi^6) / 117 +
              (245 * xi^4) / 2187 - (105 * xi^2) / 187 -
              (2450 * xi) / 6561 + 175 / 494 + fixed3(x = xi, upper = ui)
          } else if (ui == Inf | ui > xi + 1) {
            ret[i] <- fixed5(x = xi)
          } else {
            ret[i] <- 0
            }

        } else if (abs(xi == 2)) {ret[i] <- 0}
      }
      return(ret)
    },

    #' @description
    #' The squared 2-norm of the cdf is defined by
    #' \deqn{\int_a^b (F_X(u))^2 du}
    #' where X is the Distribution, \eqn{F_X} is its pdf and \eqn{a, b}
    #' are the distribution support limits.
    cdfSquared2Norm = function(x = 0, upper = 0) {


      f1 <- function(upper, x) {

        ret1 <- (1 / 26244) *
          (43379 / 33 + (28 * upper^21) / 3 - (6561 * x) / 2 -
          98 * upper^20 * x + (252315 * x^2) / 247 + (8820 * upper^19 * x^2) / 19 +
            (4900 * x^3) / 3 + (229635 * x^4) / 187 - 1701 * x^5 - 1078 * x^6 +
            (17010 * x^8) / 13 - (3605 * x^9) / 6 + (3969 * x^10) /  11 -
            (2268 * x^11) / 11 - (18 * x^15) / 143 - (7 * x^21) / 138567 -
          (280 / 3) * upper^18 * (- 1 + 14 * x^3) - 63 * upper^16 * x^2 * (- 55 + 49 * x^3) +
          (840 / 17) * upper^17 * x * (- 17 + 49 * x^3) -
            (70 / 3) * upper^3 * (- 1 + x^3)^2 * (- 280 - 729 * x^2 + 280 * x^3) +
          (630 / 13) * upper^13 * x^2 * (239 - 364 * x^3 + 14 * x^6) -
            30 * upper^14 * x * (109 - 490 * x^3 + 56 * x^6) +
          4 * upper^15 * (109 - 2170 * x^3 + 686 * x^6) +
          126 * upper^10 * x * (- 9 + 205 * x - 354 * x^4 + 30 * x^7) +
          315 * upper^4 * x * (- 27 + 70 * x + 135 * x^3 - 140 * x ^4 -
            108 * x^6 + 70 * x^7) +
          90 * upper^7 * x * (- 54 + 385 * x + 378 * x^3 - 994 * x^4 +
            105 * x^7) - (70 / 3) *
          upper^12 * (- 59 + 1101 * x ^3 - 651 * x^6 + 7 * x^9) -
          (35 / 3) * upper^9 * (- 265 - 486 * x^2 + 4580 * x^3 -
            3180 * x^6 + 80 * x^9) -
          14 * upper^6 * (- 350 - 1215 * x^2 + 4025 * x^3 + 3402 * x^5 -
            5425 * x^6 + 175 * x^9) +
          (14 / 11) * upper^11 * (162 - 6490 * x + 31155 * x^4 - 7260 * x^7 + 14 * x^10) +
          (15 / 2) * upper^8 * (162 - 1855 * x - 2268 * x^3 + 10640 * x^4 -
            3000 * x^7 + 14 * x^10) +
          21 * upper^5 * (162 - 700 * x - 1620 * x^3 + 2905 * x^4 + 2268 * x^6 -
            2300 * x^7 + 14 * x^10) +
          81 * upper * (81 - 140 * x + 105 * x^4 - 60 * x^7 + 14 * x^10) +
          70 * upper^2 * (162 - 140 * x - 243 * x^3 + 105 * x^4 + 243 * x^6 -
            60 * x^7 - 81 * x^9 + 14 * x^10))

        return(ret1)
      }


      f2 <- function(upper, x) {

        ret2 <- - ((7 * upper^21) / 19683) + (49 * upper^20 * x) / 13122 -
          (245 * upper^19 * x^2) / 13851 + (980 * upper^18 * x^3) / 19683 -
          (70 * upper^17 * x * (3 + 49 * x^3)) / 37179 +
          (7 * upper^16 * x^2 * (20 + 49 * x^3)) / 2916 +
          (35 * upper^4 * (- 1 + x) * x * (1 + x + x^2) *
             (27 - 70 * x - 108 * x^3 + 70 * x^4)) / 2916 +
          (5 * upper^14 * x * (- 11 + 350 * x^3 + 56 * x^6)) / 4374 -
          (upper^15 * (- 11 + 1190 * x^3 + 686 * x^6)) / 6561 +
          (7 * upper^10 * x * (- 9 - 5 * x + 234 * x^4 + 30 * x^7)) / 1458 -
          (7 * upper^11 * x * (- 90 + 13155 * x^3 + 7140 * x ^6 + 14 * x^9)) / 144342 -
          (35 * upper^9 * (- 55 - 486 * x^2 + 380 * x^3 + 2700 * x^6 + 80 * x ^9)) / 78732 +
          (1 / 324) * upper * (81 - 140 * x + 105 * x^4 - 60 * x^7 + 14 * x^10) -
          (35 * upper^3 * (- 1 + x^3)^2 * (- 280 + x^2 * (- 729 + 280 * x))) / 39366 -
          (35 * upper^13 * x^2 * (- 1 + 14 * x^3 * (22 + x^3))) / 18954 +
          (35 * upper^12 * x^3 * (201 + 7 * x^3 * (87 + x^3))) / 39366 +
          (5 * upper^8 * (162 - 385 * x + 2 * x^3 *
              (- 1134 + 1645 * x + 1440 * x^4 + 7 * x^7))) / 17496 -
          (7 * upper^5 * x * (420 + x^2 * (1620 - 2695 * x -
              2268 * x^3 + 2180 * x^4 + 14 * x^7))) / 8748 +
          (35 * upper^2 * (162 - 140 * x + x^3 * (- 243 + 105 * x +
              243 * x^3 - 60 * x^4 - 81 * x^6 + 14 * x^7))) / 13122 -
          (5 * upper^7 * x * (54 + 7 * x * (- 25 + x^2 *
             (- 54 + 82 * x + 15 * x^4)))) / 1458 +
          (7 * upper^6 * x^2 * (1215 + 7 * x *
            (- 425 + x^2 * (- 486 + 25 * x * (25 + x^3))))) / 13122

        return(ret2)
      }

      f3 <- function(upper, x) {

        ret3 <- (1 / 26244) *
          ((28 * upper^21) / 3 - 6561 * x - 98 * upper^20 * x +
             (8820 * upper^19 * x^2) / 19 + (9800 * x^3) / 3 +
          3402 * x^5 + (1295 * x^9) / 6 + (2268 * x^11) / 11 - (18 * x^15) / 143 -
          (7 * x^21) / 138567 - (280 / 3) * upper^18 * (1 + 14 * x^3) +
          (840 / 17) * upper^17 * x * (17 + 49 * x^3) - 63 * upper^16 * x^2 * (55 + 49 * x^3) +
          (70 / 3) * upper^3 * (1 + x^3)^2 * (280 - 729 * x^2 + 280 * x^3) +
          (630 / 13) * upper^13 * x^2 * (239 + 364 * x^3 + 14 * x^6) -
          30 * upper^14 * x * (109 + 490 * x^3 + 56 * x^6) +
          4 * upper^15 * (109 + 2170 * x^3 + 686 * x^6) -
          126 * upper^10 * x * (- 9 + 205 * x + 354 * x^4 + 30 * x^7) -
          315 * upper^4 * x * (- 27 + 70 * x - 135 * x^3 + 140 * x^4 - 108 * x^6 + 70 * x^7) +
          90 * upper^7 * x * (- 54 + 385 * x - 378 * x^3 + 994 * x^4 + 105 * x^7) -
          (70 / 3) * upper^12 * (59 + 1101 * x^3 + 651 * x^6 + 7 * x^9) +
          (35 / 3) * upper^9 * (265 - 486 * x^2 + 4580 * x^3 + 3180 * x^6 + 80 * x^9) -
          14 * upper^6 * (350 - 1215 * x^2 + 4025 * x^3 - 3402 * x^5 + 5425 * x^6 + 175 * x^9) -
          81 * upper * (- 81 + 140 * x + 105 * x^4 + 60 * x^7 + 14 * x^10) +
          21 * upper^5 * (- 162 + 700 * x - 1620 * x^3 + 2905 * x^4 - 2268 * x^6 +
              2300 * x^7 + 14 * x^10) -
          (15 / 2) * upper^8 * (- 162 + 1855 * x - 2268 * x^3 + 10640 * x^4 +
              3000 * x^7 + 14 * x^10) +
            (14 / 11) * upper^11 * (- 162 + 6490 * x + 31155 * x^4 + 7260 * x^7 + 14 * x^10) -
          70 * upper^2 * (- 162 + 140 * x - 243 * x^3 + 105 * x^4 - 243 * x^6 + 60 * x^7 -
          81 * x^9 + 14 * x^10))

        return(ret3)
      }

      f4 <- function(upper, x) {

        ret4 <- (- 140 * upper * x + 140 * x +
                   (165 * ((upper - x)^8 - (x - 1)^8) -
                    28 * ((x - 1)^11 + (upper - x)^11)) / 22 -
                    21 * ((x - 1)^5 + (upper - x)^5) + 70 * upper^2 + 81 * upper - 151) / 162

        return(ret4)
      }

      f5 <- function(upper, x) {

        ret5 <- (- ((7 * upper^21) / 19683)) + (49 * upper^20 * x) / 13122 -
          (245 * upper^19 * x^2) / 13851 +
          (980 * upper^18 * x^3) / 19683 - (70 * upper^17 * x * (3 + 49 * x^3)) / 37179 +
          (7 * upper^16 * x^2 * (20 + 49 * x^3)) / 2916 -
          (35 * upper^3 * (- 1 + x^3)^2 * (- 280 - 729 * x^2 + 280 * x^3)) / 39366 +
          (upper^15 * (11 - 1190 * x^3 - 686 * x^6)) / 6561 +
          (35 * upper^12 * x^3 * (201 + 609 * x^3 + 7 * x^6)) / 39366 -
          (35 * upper^13 * x^2 * (- 1 + 308 * x^3 + 14 * x^6)) / 18954 +
          (5 * upper^14 * x * (- 11 + 350 * x^3 + 56 * x ^6)) / 4374 +
          (7 * upper^10 * x * (- 9 - 5 * x + 234 * x^4 + 30 * x^7)) / 1458 +
          (35 * upper^4 * x * (- 27 + 70 * x + 135 * x^3 - 140 * x^4 -
                                 108 * x^6 + 70 * x^7)) / 2916 -
          (5 * upper^7 * x * (54 - 175 * x - 378 * x^3 + 574 * x^4 + 105 * x^7)) / 1458 +
          (7 * upper^6 * x^2 * (1215 - 2975 * x - 3402 * x^3 + 4375 * x^4 + 175 * x^7)) / 13122 -
          (7 * upper^5 * x * (420 + 1620 *  x^2 - 2695 * x^3 - 2268 * x^5 +
                                2180 * x^6 + 14 * x^9)) / 8748 -
          (7 * upper^11 * x * (- 90 + 13155 * x^3 + 7140 * x^6 + 14 * x^9)) / 144342 -
          (35 * upper^9 * (-55 - 486 * x^2 + 380 * x^3 + 2700 * x^6 + 80 * x^9)) / 78732 +
          (1 / 324) * upper * (81 - 140 * x + 105 * x^4 - 60 * x^7 + 14 * x^10) +
          (5 * upper^8 * (162 - 385 * x - 2268 * x^3 + 3290 * x^4 + 2880 * x^7 +
                            14 * x^10)) / 17496 +
          (35 * upper^2 * (162 - 140 * x - 243 * x^3 + 105 * x^4 + 243 * x^6 -
                             60 * x^7 - 81 * x^9 + 14 * x^10)) / 13122 +
          (3647 / 3 - (933201 * x) / 374 - (413910 * x^2) / 247 +
             6860 * x^3 - (952560 * x^4) / 187 + (45927 * x^5) / 13 -
             4802 * x^6 + (34020 * x^7) / 11 - (17010 * x^8) / 13 +
             (2135 * x^9) / 2 - (3969 * x^10) / 11 - (1085 * x^12) / 33 + (38 * x^15) / 143 -
             (35 * x^18) / 7293 + (7 * x^21) / 138567) / 26244

        return(ret5)
      }


      f6 <- function(upper, x) {

        ret6 <- (1 / 162) *
          (- 151 + 81 * upper + 70 * upper^2 +
             21 * ((upper - x)^5 + (- 1 + x)^5) +
             (15 / 2) * ((upper - x)^8 - (- 1 + x)^8) +
             (14 / 11) * ((upper - x)^11 + (- 1 + x)^11) + 140 * x - 140 * upper * x)
        return(ret6)
      }


      f7 <- function(upper, x) {
        ret7 <- (1 / 162) * (- 151 + 21 * (- 1 + x)^5 - (15 / 2) * (- 1 + x)^8 +
                               (14 / 11) * (- 1 + x)^11 + 221 * x - 70 * x^2)
        return(ret7)
      }


      f8 <- function(upper, x) {

        ret8 <- (7 * upper^21) / 19683 - (49 * upper^20 * x) / 13122 +
          (245 * upper^19 * x ^2) / 13851 - (70 * upper^18 * (- 1 + 14 * x^3)) / 19683 -
          (7 * upper^16 * x^2 * (- 55 + 49 * x^3)) / 2916 +
          (70 * upper^17 * x * (- 17 + 49 * x ^3)) / 37179 -
          (35 * upper^3 * (- 1 + x^3)^2 * (- 280 - 729 * x^2 + 280 * x^3)) / 39366 +
          (35 * upper^13 * x^2 * (239 - 364 * x^3 + 14 * x^6)) / 18954 -
          (5 * upper^14 * x * (109 - 490 * x^3 + 56 * x^6)) / 4374 +
          (upper^15 * (109 - 2170 * x^3 + 686 * x^6)) / 6561 +
          (7 * upper^10 * x * (- 9 + 205 * x - 354 * x^4 + 30 * x^7)) / 1458 +
          (35 * upper^4 * x * (- 27 + 70 * x + 135 * x^3 - 140 * x^4 -
                                 108 * x^6 + 70 * x^7)) / 2916 +
          (5 * upper^7 * x * (- 54 + 385 * x + 378 * x^3 - 994 * x^4 +
                                105 * x^7)) / 1458 -
          (35 * upper^12 * (- 59 + 1101 * x^3 - 651 * x^6 + 7 * x^9)) / 39366 -
          (35 * upper^9 * (- 265 - 486 * x^2 + 4580 * x^3 - 3180 * x^6 +
                             80 * x^9)) / 78732 -
          (7 * upper^6 * (- 350 - 1215 * x^2 + 4025 * x^3 + 3402 * x^5 -
                            5425 * x^6 + 175 * x^9)) / 13122 +
          (7 * upper^11 * (162 - 6490 * x + 31155 * x^4 - 7260 * x^7 + 14 * x^10)) / 144342 +
          (5 * upper^8 * (162 - 1855 * x - 2268 * x^3 + 10640 * x^4 -
                            3000 * x^7 + 14 * x^10)) / 17496 +
          (7 * upper^5 * (162 - 700 * x - 1620 * x^3 + 2905 * x^4 + 2268 * x^6 -
                            2300 * x^7 + 14 * x^10)) / 8748 +
          (1 / 324) * upper * (81 - 140 * x + 105 * x^4 - 60 * x^7 + 14 * x^10) +
          (35 * upper^2 * (162 - 140 * x - 243 * x^3 + 105 * x^4 + 243 * x^6 - 60 * x^7 -
                         81 * x^9 + 14 * x^10)) / 13122 +
          (43379 / 33 - (6561 * x) / 2 + (252315 * x^2) / 247 +
            (4900 * x^3) / 3 + (229635 * x^4) / 187 - 1701 * x^5 -
             1078 * x^6 + (17010 * x^8) / 13 + (2450 * x^9) / 3 +
             (3969 * x^10) / 11) / 26244

        return(ret8)
      }

      f9 <- function(upper, x) {
        ret9 <-  - ((7 * upper^21) / 19683) + (49 * upper^20 * x) / 13122 -
          (245 * upper^19 * x^2) / 13851 + (980 * upper^18 * x^3) / 19683 +
          (7 * upper^16 * x^2 * (- 20 + 49 * x^3)) / 2916 -
          (70 * upper^17 * x * (- 3 + 49 * x^3)) / 37179 +
          (35 * upper^3 * (1 + x^3)^2 * (280 - 729 * x^2 + 280 * x^3)) / 39366 +
          (upper^15 * (44 + 4760 * x^3 - 2744 * x^6)) / 26244 +
          (35 * upper^12 * x^3 * (201 - 609 * x^3 + 7 * x^6)) / 39366 -
          (35 * upper^13 * x^2 * (- 1 - 308 * x^3 + 14 * x^6)) / 18954  +
          (5 * upper^14 * x * (- 11 - 350 * x^3 + 56 * x^6)) / 4374 -
          (7 * upper^10 * x * (- 9 - 5 * x - 234 * x^4 + 30 * x^7)) / 1458 -
          (35 * upper^4 * x * (- 27 + 70 * x - 135 * x^3 + 140 * x^4 -
                                108 * x^6 + 70 * x^7)) / 2916 -
          (5 * upper^7 * x * (54 - 175 * x + 378 * x^3 - 574 * x^4 +
                                105 * x^7)) / 1458 +
          (7 * upper^6 * x^2 * (1215 - 2975 * x + 3402 * x^3 -
                                  4375 * x^4 + 175 * x^7)) / 13122 -
          (7 * upper^11 * x * (90 + 13155 * x^3 - 7140 * x^6 + 14 * x^9)) / 144342 -
          (7 * upper^5 * x * (- 420 + 1620 * x^2 - 2695 * x^3 +
                                2268 * x^5 - 2180 * x^6 + 14 * x^9)) / 8748 +
          (35 * upper^9 * (55 - 486 * x^2 + 380 * x^3 - 2700 * x^6 + 80 * x^9)) / 78732 -
          (5 * upper^8 * (- 162 + 385 * x - 2268 * x^3 + 3290 * x^4 - 2880 * x^7 +
                            14 * x^10)) / 17496 -
          (1 / 324) * upper * (- 81 + 140 * x + 105 * x^4 + 60 * x^7 + 14 * x^10) -
          (35 * upper^2 * (- 162 + 140 * x - 243 * x^3 + 105 * x^4 -
                             243 * x^6 + 60 * x^7 - 81 * x^9 + 14 * x^10)) / 13122 +
          (- 6561 * x + (9800 * x^3) / 3 + 980 * x^6 + (1505 * x^9) / 6 +
             (1085 * x^12) / 33 + (38 * x^15) / 143 + (35 * x^18) / 7293 +
             (7 * x^21) / 138567) / 26244
        return(ret9)
      }


      f10 <- function(upper, x) {
        ret10 <- - ((7 * upper^21) / 19683) + (49 * upper^20 * x) / 13122 -
          (245 * upper^19 * x^2) / 13851 + (70 * upper^18 * (1 + 14 * x^3)) / 19683 -
          (70 * upper^17 * x * (17 + 49 * x^3)) / 37179 +
          (7 * upper^16 * x^2 * (55 + 49 * x^3)) / 2916 -
          (35 * upper^4 * x * (1 + x) * (1 + (- 1 + x) * x) *
             (- 27 + 70 * x - 108 * x^3 + 70 * x^4)) / 2916 +
          (5 * upper^14 * x * (39 + 490 * x^3 + 56 * x^6)) / 4374 -
          (2 * upper^15 * (30 + 1085 * x^3 + 343 * x^6)) / 6561 +
          (7 * upper^10 * x * (9 - 145 * x - 234 * x^4 + 30 * x^7))  / 1458 -
          (7 * upper^11 * x * (- 2610 - 12945 * x^3 + 7260 * x^6 + 14 * x^9)) / 144342 -
          (35 * upper^9 * (- 105 + 486 * x^2 - 4100 * x^3 - 2700 * x^6 + 80 * x^9)) / 78732 -
          (1 / 324) * upper * (- 81 + 140 * x + 105 * x^4 + 60 * x^7 + 14 * x^10) +
          (35 * upper^3 * (1 + x^3)^2 * (280 + x^2 * (- 729 + 280 * x))) / 39366 -
          (35 * upper^13 * x^2 * (29 + 14 * x^3 * (26 + x^3))) / 18954 +
          (35 * upper^12 * x^3 * (- 159 + 7 * x^3 * (93 + x^3))) / 39366 +
          (5 * upper^8 * x * (- 1575 + 2 * x^2 * (1134 - 5215 * x - 1440 * x^4 +
                7 * x^7))) / 17496 +
          (7 * upper^5 * (- 162 + 700 * x + x^3 * (- 1620 + 2905 * x -
                2268 * x^3 + 2300 * x^4 + 14 * x^7))) / 8748 -
          (35 * upper^2 * (- 162 + 140 * x + x^3 * (- 243 + 105 * x -
                243 * x^3 + 60 * x^4 - 81 * x^6 + 14 * x^7))) / 13122 +
          (5 * upper^7 * x * (-54 + 7 * x * (55 + x^2 * (- 54 + 142 * x +
                15 * x^4)))) / 1458 -
          (7 * upper^6 * (350 + x^2 * (- 1215 + 7 * x * (575 + x^2 *
                (- 486 + 25 * x * (31 + x^3)))))) / 13122

        return(ret10)
      }

      f11 <- function(upper, x) {
        ret11 <- (- 3322 + 1782 * upper + 1540 * upper^2 - 4862 * x - 1540 * x^2 -
                 462 * (upper^5 - (1 + x)^5) - 165 * (upper^8 - (1 + x)^8) +
                 28 * (upper^11 - (1 + x)^11)) / 3564
        return(ret11)
      }

      f12 <- function(upper, x) {
        ret12 <- - ((7 * upper^21) / 19683) + (49 * upper^20 * x) / 13122 -
          (245 * upper^19 * x^2) / 13851 +
          (980 * upper^18 * x^3) / 19683 + (7 * upper^16 * x^2 * (- 20 + 49 * x^3)) / 2916 -
          (70 * upper^17 * x * (- 3 + 49 * x^3)) / 37179 +
          (35 * upper^3 * (1 + x^3)^2 * (280 - 729 * x^2 + 280 * x^3)) / 39366 +
          (upper^15 * (44 + 4760 * x^3 - 2744 * x^6)) / 26244 +
          (35 * upper^12 * x^3 * (201 - 609 * x^3 + 7 * x^6)) / 39366 -
          (35 * upper^13 * x^2 * (- 1 - 308 * x^3 + 14 * x^6)) / 18954 +
          (5 * upper^14 * x * (- 11 - 350 * x^3 + 56 * x^6)) / 4374 -
          (7 * upper^10 * x * (- 9 - 5 * x - 234 * x^4 + 30 * x^7)) / 1458 -
          (35 * upper^4 * x * (- 27 + 70 * x - 135 * x^3 + 140 * x^4 -
            108 * x^6 + 70 * x^7)) / 2916 -
          (5 * upper^7 * x * (54 - 175 * x + 378 * x^3 - 574 * x^4 + 105 * x^7)) / 1458 +
          (7 * upper^6 * x^2 * (1215 - 2975 * x + 3402 * x^3 - 4375 * x^4 + 175 * x^7)) / 13122 -
          (7 * upper^11 * x * (90 + 13155 * x^3 - 7140 * x^6 + 14 * x^9)) / 144342 -
          (7 * upper^5 * x * (- 420 + 1620 * x^2 - 2695 * x^3 + 2268 * x^5 -
            2180 * x^6 + 14 * x^9)) / 8748 +
          (35 * upper^9 * (55 - 486 * x^2 + 380 * x^3 - 2700 * x^6 + 80 * x^9)) / 78732 +
          (3647 / 3 - (1520613 * x) / 374 - (413910 * x^2) / 247 - (10780 * x^3) / 3 -
             (952560 * x^4) / 187 - (45927 * x^5) / 13 -
             3822 * x^6 - (34020 * x^7) / 11 - (17010 * x^8) / 13 - (2450 * x^9) / 3 -
             (3969 * x^10) / 11) / 26244 -
          (5 * upper^8 * (- 162 + 385 * x - 2268 * x^3 + 3290 * x^4 - 2880 * x^7 +
            14 * x^10)) / 17496 -
          (1 / 324) * upper * (- 81 + 140 * x + 105 * x^4 + 60 * x^7 + 14 * x^10) -
          (35 * upper^2 * (- 162 + 140 * x - 243 * x^3 + 105 * x^4 - 243 * x^6 + 60 * x^7 -
                         81 * x^9 + 14 * x^10)) / 13122

        return(ret12)
      }


      f13 <- function(upper, x) {

        ret13 <- (- 3322 + 1782 * upper + 1540 * upper^2 - 4862 * x - 1540 * x ^2 +
                 462 * (upper^5 - (1 + x)^5) + 165 * (upper^8 - (1 + x)^8) +
                 28 * (upper^11 - (1 + x)^11)) / 3564
        return(ret13)
      }


      f14 <- function(upper, x) {
        ret14 <- upper / 2 + (35 * upper^2) / 81 - (7 * upper^5) / 54 +
          (5 * upper^8) / 108 - (7 * upper^11) / 891
        return(ret14)
      }


      xl <- length(x)
      ul <- length(upper)
      len <- max(xl, ul)

      ret <- numeric(len)
      for (i in seq(len)) {

        xi <- x[ifelse(i %% xl == 0, xl, i %% xl)]
        ui <- upper[ifelse(i %% ul == 0, ul, i %% ul)]

        if (xi >= 0 & xi <= 1) {
          if (ui <= xi - 1) {
            ret[i] <- 0
          } else if (ui >= xi - 1 & ui <= 0) {
            ret[i] <- f1(upper = ui, x = xi)
          } else if (ui >= 0 & ui <= xi) {
            ret[i] <- f1(upper  = 0, x = xi) + f2(upper = ui, x = xi)
          } else if (ui >= xi & ui <= 1) {
            ret[i] <- f1(upper  = 0, x = xi) + f2(upper = xi, x = xi) +
              f3(upper = ui, x  = xi)
          } else if (ui >= 1 & ui <= xi + 1) {
            ret[i] <- f1(upper  = 0, x = xi) + f2(upper = xi, x = xi) +
              f3(upper = 1, x = xi) + f4(upper = ui, x = xi)
          } else if (ui >= xi + 1) {
            ret[i] <- f1(upper  = 0, x = xi) + f2(upper = xi, x = xi) +
              f3(upper = 1, x = xi) + f4(upper = xi + 1, x = xi) +
              (ui - xi - 1)
          }
        } else if (xi >= 1 & xi <= 2) {
          if (ui <= xi - 1) {
            ret[i] <- 0
          } else if (ui >= xi - 1 & ui <= 1) {
            ret[i] <- f5(upper = ui, x = xi)
          } else if (ui >= 1 & ui <= xi) {
            ret[i] <- f5(upper = 1, x = xi) + f6(upper = ui, x = xi)
          } else if (ui >= xi & ui <= xi + 1) {
            ret[i] <- f5(upper = 1, x = xi) + f6(upper = xi, x = xi) +
              f7(upper = ui, x = xi)
          } else if (ui >= xi + 1) {
            ret[i] <- f5(upper = 1, x = xi) + f6(upper = xi, x = xi) +
              f7(upper = xi + 1, x = xi) + (ui - xi - 1)
          }
        } else if (xi >= 2) {
          if (ui <= xi - 1) {
            ret[i] <- 0
          } else if (ui >= xi - 1 & ui <= xi) {
            ret[i] <- (1 / 162) *
              (567 / 22 + 81 * ui + 70 * ui^2 + 21 * (ui - xi)^5 +
                 (15 / 2) * (ui - xi)^8 + (14 / 11) * (ui - xi)^11 -
                 81 * xi - 140 * ui * xi + 70 * xi^2)
          } else if (ui >= xi & ui <= xi + 1) {
            ret[i] <- 7 / 44 + (1 / 162) *
              (81 * ui + 70 * ui^2 - 21 * (ui - xi)^5 +
                 (15 / 2) * (ui - xi)^8 - (14 / 11) * (ui - xi)^11 -
                 81 * xi - 140 * ui * xi + 70 * xi^2)
          } else if (ui >= xi + 1) {
            ret[i] <- ui - xi
          }
        } else if (xi >= - 1 & xi <= 0) {
          if (ui <= - 1) {
            ret[i] <- 0
          } else if (ui >= - 1 & ui <= xi) {
            ret[i] <- f8(upper = ui, x = xi)
          } else if (ui >= xi & ui <= 0) {
            ret[i] <- f8(upper = xi, x = xi) + f9(upper = xi, x = xi)
          } else if (ui >= 0 & ui <= xi + 1) {
            ret[i] <- f8(upper = xi, x = xi) + f9(upper = 0, x = xi) +
              f10(upper = ui, x = xi)
          } else if (ui >= xi + 1 & ui <= 1) {
            ret[i] <- f8(upper = xi, x = xi) + f9(upper = 0, x = xi) +
              f10(upper = xi + 1, x = xi) + f11(upper = ui, x = xi)
          } else if (ui >= 1) {
            ret[i] <-  f8(upper = xi, x = xi) + f9(upper = 0, x = xi) +
              f10(upper = xi + 1, x = xi) + f11(upper = 1, x = xi) + (ui - 1)
          }
        } else if (xi >= - 2 & xi <= - 1) {
          if (ui <= - 1) {
            ret[i] <- 0
          } else if (ui >= - 1 & ui <= xi + 1) {
            ret[i] <- f12(upper = ui, x = xi)
          } else if (ui >= xi + 1 & ui <= 0) {
            ret[i] <- f12(upper = xi + 1, x = xi) + f13(upper = ui, x = xi)
          } else if (ui >= 0 & ui <= 1) {
            ret[i] <- f12(upper = xi + 1, x = xi) + f13(upper = 0, x = xi) +
              f14(upper = ui, x = xi)
          } else if (ui >= 1) {
            ret[i] <- f12(upper = xi + 1, x = xi) + f13(upper = 0, x = xi) +
              f14(upper = 1, x = xi)
          }
        } else if (xi <= - 2) {
          if (ui <= - 1) {
            ret[i] <- 0
          } else if (ui >= - 1 & ui <= 0) {
            ret[i] <- 7 / 44 + ui / 2 + (35 * ui^2) / 81 +
              (7 * ui^5) / 54 + (5 * ui^8) / 108 +
              (7 * ui^11) / 891
          } else if (ui >= 0 & ui <= 1) {
            ret[i] <- 7 / 44 + ui / 2 + (35 * ui^2) / 81 -
              (7 * ui^5) / 54 + (5 * ui^8) / 108 - (7 * ui^11) / 891
          } else if (ui >= 1) {
            ret[i] <- ui
          }
        }
      }
      return(ret)
    },

    #' @description
    #' The variance of a distribution is defined by the formula
    #' \deqn{var_X = E[X^2] - E[X]^2}
    #' where \eqn{E_X} is the expectation of distribution X. If the distribution is multivariate the
    #' covariance matrix is returned.
    #' @param ... Unused.
    variance = function(...) {
      return(35 / 243)
    }
  ),

  private = list(
    .isQuantile = 0L,
    .pdf = function(x, log = FALSE) {
      C_TricubeKernelPdf(x, log)
    },
    .cdf = function(x, lower.tail = TRUE, log.p = FALSE) {
      C_TricubeKernelCdf(x, lower.tail, log.p)
    }
  )
)

.distr6$kernels <- rbind(
  .distr6$kernels,
  data.table::data.table(
    ShortName = "Tric", ClassName = "Tricube",
    Support = "[-1,1]", Packages = "-"
  )
)
RaphaelS1/distr6 documentation built on Feb. 24, 2024, 9:14 p.m.