R/poly_skurt_check.R

In SimMultiCorrData: Simulation of Correlated Data with Multiple Variable Types

#' @title Headrick's Fifth-Order Transformation Lagrangean Constraints for Lower Boundary of Standardized Kurtosis
#'
#' @description This function gives the first-order conditions of the multi-constraint Lagrangean expression
#'     \deqn{F(c1, ..., c5, \lambda_{1}, ..., \lambda_{4}) = f(c1, ..., c5) + \lambda_{1} * [1 - g(c1, ..., c5)]}
#'     \deqn{+ \lambda_{2} * [\gamma_{1} - h(c1, ..., c5)] + \lambda_{3} * [\gamma_{3} - i(c1, ..., c5)]}
#'     \deqn{+ \lambda_{4} * [\gamma_{4} - j(c1, ..., c5)]}
#'     used to find the lower kurtosis boundary for a given skewness and standardized fifth and sixth cumulants
#'     in \code{\link[SimMultiCorrData]{calc_lower_skurt}}.  The partial derivatives are described in Headrick (2002,
#'     \doi{10.1016/S0167-9473(02)00072-5}), but he does not provide
#'     the actual equations. The equations used here were found with \code{D} (see \code{\link[stats]{deriv}}).
#'     Here, \eqn{\lambda_{1}, ..., \lambda_{4}} are the Lagrangean multipliers, \eqn{\gamma_{1}, \gamma_{3}, \gamma_{4}} are the user-specified
#'     values of skewness, fifth cumulant, and sixth cumulant, and \eqn{f, g, h, i, j} are the equations for standardized kurtosis, variance,
#'     fifth cumulant, and sixth cumulant expressed in terms of the constants.  This function would not ordinarily be called by the user.
#'
#' @param c a vector of constants c1, ..., c5, lambda1, ..., lambda4
#' @param a a vector of skew, fifth standardized cumulant, sixth standardized cumulant
#' @export
#' @keywords kurtosis, boundary, Headrick
#' @seealso \code{\link[SimMultiCorrData]{calc_lower_skurt}}
#' @return A list with components:
#' @return \eqn{dF/d\lambda_{1} = 1 - g(c1, ..., c5)}
#' @return \eqn{dF/d\lambda_{2} = \gamma_{1} - h(c1, ..., c5)}
#' @return \eqn{dF/d\lambda_{3} = \gamma_{3} - i(c1, ..., c5)}
#' @return \eqn{dF/d\lambda_{4} = \gamma_{4} - j(c1, ..., c5)}
#' @return \eqn{dF/dc1 = df/dc1 - \lambda_{1} * dg/dc1 - \lambda_{2} * dh/dc1 - \lambda_{3} * di/dc1 - \lambda_{4} * dj/dc1}
#' @return \eqn{dF/dc2 = df/dc2 - \lambda_{1} * dg/dc2 - \lambda_{2} * dh/dc2 - \lambda_{3} * di/dc2 - \lambda_{4} * dj/dc2}
#' @return \eqn{dF/dc3 = df/dc3 - \lambda_{1} * dg/dc3 - \lambda_{2} * dh/dc3 - \lambda_{3} * di/dc3 - \lambda_{4} * dj/dc3}
#' @return \eqn{dF/dc4 = df/dc4 - \lambda_{1} * dg/dc4 - \lambda_{2} * dh/dc4 - \lambda_{3} * di/dc4 - \lambda_{4} * dj/dc4}
#' @return \eqn{dF/dc5 = df/dc5 - \lambda_{1} * dg/dc5 - \lambda_{2} * dh/dc5 - \lambda_{3} * di/dc5 - \lambda_{4} * dj/dc5}
#' @return If the suppled values for \code{c} and \code{a} satisfy the Lagrangean expression, it will return 0 for each component.
#' @references
#' Headrick TC (2002). Fast Fifth-order Polynomial Transforms for Generating Univariate and Multivariate
#'     Non-normal Distributions. Computational Statistics & Data Analysis, 40(4):685-711. \doi{10.1016/S0167-9473(02)00072-5}.
#'     (\href{http://www.sciencedirect.com/science/article/pii/S0167947302000725}{ScienceDirect})
#'
#' Headrick TC (2004). On Polynomial Transformations for Simulating Multivariate Nonnormal Distributions.
#'     Journal of Modern Applied Statistical Methods, 3(1), 65-71. \doi{10.22237/jmasm/1083370080}.
#'
#' Headrick TC, Kowalchuk RK (2007). The Power Method Transformation: Its Probability Density Function, Distribution
#'     Function, and Its Further Use for Fitting Data. Journal of Statistical Computation and Simulation, 77, 229-249. \doi{10.1080/10629360600605065}.
#'
#' Headrick TC, Sheng Y, & Hodis FA (2007). Numerical Computing and Graphics for the Power Method Transformation Using
#'     Mathematica. Journal of Statistical Software, 19(3), 1 - 17. \doi{10.18637/jss.v019.i03}.
#'
poly_skurt_check <- function(c, a) {
  F <- numeric(9)
  g1 <- a[1]
  g3 <- a[2]
  g4 <- a[3]
  f <- 24 * (2 * c[2]^4 + 96 * c[2]^3 * c[4] + c[1]^3 * (c[3] + 10 * c[5]) +
             30 * c[2]^2 * (6 * c[3]^2 + 64 * c[4]^2 + 140 * c[3] * c[5] +
             945 * c[5]^2) + c[1]^2 * (2 * c[2]^2 + 18 * c[3]^2 + 36 * c[2] *
             c[4] + 192 * c[4]^2 + 375 * c[3] * c[5] + 2250 * c[5]^2) + 36 *
             c[2] * c[4] * (125 * c[3]^2 + 528 * c[4]^2 + 3360 * c[3] * c[5] +
             25725 * c[5]^2) + 3 * c[1] * (45 * c[3]^3 + 1584 * c[3] *
             c[4]^2 + 1590 * c[3]^2 * c[5] + 21360 * c[4]^2 * c[5] + 21525 *
             c[3] * c[5]^2 + 110250 * c[5]^3 + 12 * c[2]^2 * (c[3] +
             10 * c[5]) + 8 * c[2] * c[4] * (32 * c[3] + 375 * c[5])) + 9 *
             (45 * c[3]^4 + 8704 * c[4]^4 + 2415 * c[3]^3 * c[5] + 932400 *
             c[4]^2 * c[5]^2 + 3018750 * c[5]^4 + 20 * c[3]^2 * (178 *
             c[4]^2 + 2765 * c[5]^2) + 35 * c[3] * (3104 * c[4]^2 * c[5] +
                                                      18075 * c[5]^3)))
  f_c1 <- 24 * (3 * c[1]^2 * (c[3] + 10 * c[5]) + 2 * c[1] * (2 * c[2]^2 +
                18 * c[3]^2 + 36 * c[2] * c[4] + 192 * c[4]^2 + 375 * c[3] *
                c[5] + 2250 * c[5]^2) + 3 * (45 * c[3]^3 + 1584 * c[3] *
                c[4]^2 + 1590 * c[3]^2 * c[5] + 21360 * c[4]^2 * c[5] +
                21525 * c[3] * c[5]^2 + 110250 * c[5]^3 + 12 * c[2]^2 *
                (c[3] + 10 * c[5]) + 8 * c[2] * c[4] * (32 * c[3] + 375 *
                                                          c[5])))
  f_c2 <- 24 * (2 * (4 * c[2]^3) + 96 * (3 * c[2]^2) * c[4] + 30 * (2 *
               c[2]) * (6 * c[3]^2 + 64 * c[4]^2 + 140 * c[3] * c[5] +
               945 * c[5]^2) + c[1]^2 * (2 * (2 * c[2]) + 36 * c[4]) + 36 *
               c[4] * (125 * c[3]^2 + 528 * c[4]^2 + 3360 * c[3] * c[5] +
               25725 * c[5]^2) + 3 * c[1] * (12 * (2 * c[2]) * (c[3] + 10 *
               c[5]) + 8 * c[4] * (32 * c[3] + 375 * c[5])))
  f_c3 <- 24 * (c[1]^3 + 30 * c[2]^2 * (6 * (2 * c[3]) + 140 * c[5]) +
                c[1]^2 * (18 * (2 * c[3]) + 375 * c[5]) + 36 * c[2] * c[4] *
                (125 * (2 * c[3]) + 3360 * c[5]) + 3 * c[1] * (45 * (3 *
                c[3]^2) + 1584 * c[4]^2 + 1590 * (2 * c[3]) * c[5] + 21525 *
                c[5]^2 + 12 * c[2]^2 + 8 * c[2] * c[4] * 32) + 9 * (45 * (4 *
                c[3]^3) + 2415 * (3 * c[3]^2) * c[5] + 20 * (2 * c[3]) *
                (178 * c[4]^2 + 2765 * c[5]^2) + 35 * (3104 * c[4]^2 * c[5] +
                                                         18075 * c[5]^3)))
  f_c4 <- 24 * (96 * c[2]^3 + 30 * c[2]^2 * (64 * (2 * c[4])) + c[1]^2 * (36 *
                c[2] + 192 * (2 * c[4])) + (36 * c[2] * (125 * c[3]^2 + 528 *
                c[4]^2 + 3360 * c[3] * c[5] + 25725 * c[5]^2) + 36 * c[2] *
                c[4] * (528 * (2 * c[4]))) + 3 * c[1] * (1584 * c[3] * (2 *
                c[4]) + 21360 * (2 * c[4]) * c[5] + 8 * c[2] * (32 * c[3] +
                375 * c[5])) + 9 * (8704 * (4 * c[4]^3) + 932400 * (2 *
                c[4]) * c[5]^2 + 20 * c[3]^2 * (178 * (2 * c[4])) + 35 *
                c[3] * (3104 * (2 * c[4]) * c[5])))
  f_c5 <- 24 * (c[1]^3 * 10 + 30 * c[2]^2 * (140 * c[3] + 945 * (2 * c[5])) +
                c[1]^2 * (375 * c[3] + 2250 * (2 * c[5])) + 36 * c[2] * c[4] *
                (3360 * c[3] + 25725 * (2 * c[5])) + 3 * c[1] * (1590 *
                c[3]^2 + 21360 * c[4]^2 + 21525 * c[3] * (2 * c[5]) + 110250 *
                (3 * c[5]^2) + 12 * c[2]^2 * 10 + 8 * c[2] * c[4] * 375) + 9 *
                (2415 * c[3]^3 + 932400 * c[4]^2 * (2 * c[5]) + 3018750 * (4 *
                 c[5]^3) + 20 * c[3]^2 * (2765 * (2 * c[5])) + 35 * c[3] *
                (3104 * c[4]^2 + 18075 * (3 * c[5]^2))))
  g <- c[1]^2 + 2 * c[2]^2 + 24 * c[2] * c[4] + 6 * c[1] * (c[3] + 5 * c[5]) +
    3 * (5 * c[3]^2 + 32 * c[4]^2 + 70 * c[3] * c[5] + 315 * c[5]^2)
  g_c1 <- 2 * c[1] + 6 * (c[3] + 5 * c[5])
  g_c2 <- 2 * (2 * c[2]) + 24 * c[4]
  g_c3 <- 6 * c[1] + 3 * (5 * (2 * c[3]) + 70 * c[5])
  g_c4 <- 24 * c[2] + 3 * (32 * (2 * c[4]))
  g_c5 <- 6 * c[1] * 5 + 3 * (70 * c[3] + 315 * (2 * c[5]))
  h <- 2 * (4 * c[2]^3 + 108 * c[2]^2 * c[4] + 3 * c[1]^2 * (c[2] + 6 *
            c[4]) + 18 * c[1] * (2 * c[2] * c[3] + 16 * c[3] * c[4] + 15 *
            c[2] * c[5] + 150 * c[4] * c[5]) + 9 * c[2] * (15 * c[3]^2 +
            128 * c[4]^2 + 280 * c[3] * c[5] + 1575 * c[5]^2) + 54 * c[4] *
            (25 * c[3]^2 + 88 * c[4]^2 + 560 * c[3] * c[5] + 3675 * c[5]^2))
  h_c1 <- 2 * (3 * (2 * c[1]) * (c[2] + 6 * c[4]) + 18 * (2 * c[2] * c[3] +
              16 * c[3] * c[4] + 15 * c[2] * c[5] + 150 * c[4] * c[5]))
  h_c2 <- 2 * (4 * (3 * c[2]^2) + 108 * (2 * c[2]) * c[4] + 3 * c[1]^2 +
              18 * c[1] * (2 * c[3] + 15 * c[5]) + 9 * (15 * c[3]^2 + 128 *
              c[4]^2 + 280 * c[3] * c[5] + 1575 * c[5]^2))
  h_c3 <- 2 * (18 * c[1] * (2 * c[2] + 16 * c[4]) + 9 * c[2] * (15 * (2 *
              c[3]) + 280 * c[5]) + 54 * c[4] * (25 * (2 * c[3]) + 560 * c[5]))
  h_c4 <- 2 * (108 * c[2]^2 + 3 * c[1]^2 * 6 + 18 * c[1] * (16 * c[3] + 150 *
              c[5]) + 9 * c[2] * (128 * (2 * c[4])) + (54 * (25 * c[3]^2 +
              88 * c[4]^2 + 560 * c[3] * c[5] + 3675 * c[5]^2) + 54 *
                c[4] * (88 * (2 * c[4]))))
  h_c5 <- 2 * (18 * c[1] * (15 * c[2] + 150 * c[4]) + 9 * c[2] * (280 * c[3] +
            1575 * (2 * c[5])) + 54 * c[4] * (560 * c[3] + 3675 * (2 * c[5])))
  i <- 24 * (16 * c[2]^5 + 5 * c[1]^4 * c[4] + 1200 * c[2]^4 * c[4] + 10 *
             c[1]^3 * (3 * c[2] * c[3] + 42 * c[3] * c[4] + 40 * c[2] * c[5] +
             570 * c[4] * c[5]) + 300 * c[2]^3 * (10 * c[3]^2 + 128 * c[4]^2 +
             280 * c[3] * c[5] + 2205 * c[5]^2) + 1080 * c[2]^2 * c[4] *
             (125 * c[3]^2 + 3920 * c[3] * c[5] + 28 * (22 * c[4]^2 + 1225 *
              c[5]^2)) + 10 * c[1]^2 * (2 * c[2]^3 + 72 * c[2]^2 * c[4] + 3 *
              c[2] * (24 * c[3]^2 + 320 * c[4]^2 + 625 * c[3] * c[5] + 4500 *
              c[5]^2) + 9 * c[4] * (109 * c[3]^2 + 528 * c[4]^2 + 3130 *
              c[3] * c[5] + 24975 * c[5]^2)) + 30 * c[1] * (8 * c[2]^3 * (2 *
              c[3] + 25 * c[5]) + 40 * c[2]^2 * c[4] * (16 * c[3] + 225 *
              c[5]) + 3 * c[2] * (75 * c[3]^3 + 3168 * c[3] * c[4]^2 + 3180 *
              c[3]^2 * c[5] + 49840 * c[4]^2 * c[5] + 50225 * c[3] * c[5]^2 +
              294000 * c[5]^3) + 6 * c[4] * (555 * c[3]^3 + 8704 * c[3] *
              c[4]^2 + 26225 * c[3]^2 * c[5] + 152160 * c[4]^2 * c[5] +
              459375 * c[3] * c[5]^2 + 2963625 * c[5]^3)) + 90 * c[2] * (270 *
              c[3]^4 + 16905 * c[3]^3 * c[5] + 280 * c[3]^2 * (89 * c[4]^2 +
              1580 * c[5]^2) + 35 * c[3] * (24832 * c[4]^2 * c[5] + 162675 *
              c[5]^3) + 4 * (17408 * c[4]^4 + 2097900 * c[4]^2 * c[5]^2 +
              7546875 * c[5]^4)) + 27 * c[4] * (14775 * c[3]^4 + 1028300 *
              c[3]^3 * c[5] + 50 * c[3]^2 * (10144 * c[4]^2 + 594055 *
              c[5]^2) + 700 * c[3] * (27904 * c[4]^2 * c[5] + 598575 *
              c[5]^3) + 3 * (316928 * c[4]^4 + 68908000 * c[4]^2 * c[5]^2 +
                               806378125 * c[5]^4)))
  i_c1 <- 24 * (5 * (4 * c[1]^3) * c[4] + 10 * (3 * c[1]^2) * (3 * c[2] *
                c[3] + 42 * c[3] * c[4] + 40 * c[2] * c[5] + 570 * c[4] *
                c[5]) + 10 * (2 * c[1]) * (2 * c[2]^3 + 72 * c[2]^2 * c[4] +
                3 * c[2] * (24 * c[3]^2 + 320 * c[4]^2 + 625 * c[3] * c[5] +
                4500 * c[5]^2) + 9 * c[4] * (109 * c[3]^2 + 528 * c[4]^2 +
                3130 * c[3] * c[5] + 24975 * c[5]^2)) + 30 * (8 * c[2]^3 *
                (2 * c[3] + 25 * c[5]) + 40 * c[2]^2 * c[4] * (16 * c[3] +
                225 * c[5]) + 3 * c[2] * (75 * c[3]^3 + 3168 * c[3] * c[4]^2 +
                3180 * c[3]^2 * c[5] + 49840 * c[4]^2 * c[5] + 50225 * c[3] *
                c[5]^2 + 294000 * c[5]^3) + 6 * c[4] * (555 * c[3]^3 + 8704 *
                c[3] * c[4]^2 + 26225 * c[3]^2 * c[5] + 152160 * c[4]^2 *
                  c[5] + 459375 * c[3] * c[5]^2 + 2963625 * c[5]^3)))
  i_c2 <- 24 * (16 * (5 * c[2]^4) + 1200 * (4 * c[2]^3) * c[4] + 10 * c[1]^3 *
                (3 * c[3] + 40 * c[5]) + 300 * (3 * c[2]^2) * (10 * c[3]^2 +
                128 * c[4]^2 + 280 * c[3] * c[5] + 2205 * c[5]^2) + 1080 *
                (2 * c[2]) * c[4] * (125 * c[3]^2 + 3920 * c[3] * c[5] + 28 *
                (22 * c[4]^2 + 1225 * c[5]^2)) + 10 * c[1]^2 * (2 * (3 *
                c[2]^2) + 72 * (2 * c[2]) * c[4] + 3 * (24 * c[3]^2 + 320 *
                c[4]^2 + 625 * c[3] * c[5] + 4500 * c[5]^2)) + 30 * c[1] *
                (8 * (3 * c[2]^2) * (2 * c[3] + 25 * c[5]) + 40 * (2 * c[2]) *
                c[4] * (16 * c[3] + 225 * c[5]) + 3 * (75 * c[3]^3 + 3168 *
                c[3] * c[4]^2 + 3180 * c[3]^2 * c[5] + 49840 * c[4]^2 * c[5] +
                50225 * c[3] * c[5]^2 + 294000 * c[5]^3)) + 90 * (270 *
                c[3]^4 + 16905 * c[3]^3 * c[5] + 280 * c[3]^2 * (89 * c[4]^2 +
                1580 * c[5]^2) + 35 * c[3] * (24832 * c[4]^2 * c[5] + 162675 *
                c[5]^3) + 4 * (17408 * c[4]^4 + 2097900 * c[4]^2 * c[5]^2 +
                                 7546875 * c[5]^4)))
  i_c3 <- 24 * (10 * c[1]^3 * (3 * c[2] + 42 * c[4]) + 300 * c[2]^3 * (10 *
                (2 * c[3]) + 280 * c[5]) + 1080 * c[2]^2 * c[4] * (125 *
                (2 * c[3]) + 3920 * c[5]) + 10 * c[1]^2 * (3 * c[2] * (24 *
                (2 * c[3]) + 625 * c[5]) + 9 * c[4] * (109 * (2 * c[3]) +
                 3130 * c[5])) + 30 * c[1] * (8 * c[2]^3 * 2 + 40 * c[2]^2 *
                c[4] * 16 + 3 * c[2] * (75 * (3 * c[3]^2) + 3168 * c[4]^2 +
                3180 * (2 * c[3]) * c[5] + 50225 * c[5]^2) + 6 * c[4] * (555 *
                (3 * c[3]^2) + 8704 * c[4]^2 + 26225 * (2 * c[3]) * c[5] +
                  459375 * c[5]^2)) + 90 * c[2] * (270 * (4 * c[3]^3) +
                16905 * (3 * c[3]^2) * c[5] + 280 * (2 * c[3]) * (89 *
                c[4]^2 + 1580 * c[5]^2) + 35 * (24832 * c[4]^2 * c[5] +
                162675 * c[5]^3)) + 27 * c[4] * (14775 * (4 * c[3]^3) +
                1028300 * (3 * c[3]^2) * c[5] + 50 * (2 * c[3]) * (10144 *
                c[4]^2 + 594055 * c[5]^2) + 700 * (27904 * c[4]^2 * c[5] +
                                                     598575 * c[5]^3)))
  i_c4 <- 24 * (5 * c[1]^4 + 1200 * c[2]^4 + 10 * c[1]^3 * (42 * c[3] + 570 *
                c[5]) + 300 * c[2]^3 * (128 * (2 * c[4])) + (1080 * c[2]^2 *
               (125 * c[3]^2 + 3920 * c[3] * c[5] + 28 * (22 * c[4]^2 + 1225 *
                c[5]^2)) + 1080 * c[2]^2 * c[4] * (28 * (22 * (2 * c[4])))) +
                10 * c[1]^2 * (72 * c[2]^2 + 3 * c[2] * (320 * (2 * c[4])) +
                (9 * (109 * c[3]^2 + 528 * c[4]^2 + 3130 * c[3] * c[5] +
                24975 * c[5]^2) + 9 * c[4] * (528 * (2 * c[4])))) + 30 *
                 c[1] * (40 * c[2]^2 * (16 * c[3] + 225 * c[5]) + 3 * c[2] *
                (3168 * c[3] * (2 * c[4]) + 49840 * (2 * c[4]) * c[5]) + (6 *
                (555 * c[3]^3 + 8704 * c[3] * c[4]^2 + 26225 * c[3]^2 * c[5] +
                152160 * c[4]^2 * c[5] + 459375 * c[3] * c[5]^2 + 2963625 *
                c[5]^3) + 6 * c[4] * (8704 * c[3] * (2 * c[4]) + 152160 * (2 *
                c[4]) * c[5]))) + 90 * c[2] * (280 * c[3]^2 * (89 * (2 *
                c[4])) + 35 * c[3] * (24832 * (2 * c[4]) * c[5]) + 4 *
                (17408 * (4 * c[4]^3) + 2097900 * (2 * c[4]) * c[5]^2)) +
                (27 * (14775 * c[3]^4 + 1028300 * c[3]^3 * c[5] + 50 *
                c[3]^2 * (10144 * c[4]^2 + 594055 * c[5]^2) + 700 * c[3] *
                (27904 * c[4]^2 * c[5] + 598575 * c[5]^3) + 3 * (316928 *
                c[4]^4 + 68908000 * c[4]^2 * c[5]^2 + 806378125 * c[5]^4)) +
                27 * c[4] * (50 * c[3]^2 * (10144 * (2 * c[4])) + 700 *
                c[3] * (27904 * (2 * c[4]) * c[5]) + 3 * (316928 * (4 *
                          c[4]^3) + 68908000 * (2 * c[4]) * c[5]^2))))
  i_c5 <- 24 * (10 * c[1]^3 * (40 * c[2] + 570 * c[4]) + 300 * c[2]^3 * (280 *
                c[3] + 2205 * (2 * c[5])) + 1080 * c[2]^2 * c[4] * (3920 *
                c[3] + 28 * (1225 * (2 * c[5]))) + 10 * c[1]^2 * (3 * c[2] *
               (625 * c[3] + 4500 * (2 * c[5])) + 9 * c[4] * (3130 * c[3] +
                24975 * (2 * c[5]))) + 30 * c[1] * (8 * c[2]^3 * 25 + 40 *
                c[2]^2 * c[4] * 225 + 3 * c[2] * (3180 * c[3]^2 + 49840 *
                c[4]^2 + 50225 * c[3] * (2 * c[5]) + 294000 * (3 * c[5]^2)) +
                6 * c[4] * (26225 * c[3]^2 + 152160 * c[4]^2 + 459375 * c[3] *
                (2 * c[5]) + 2963625 * (3 * c[5]^2))) + 90 * c[2] * (16905 *
                c[3]^3 + 280 * c[3]^2 * (1580 * (2 * c[5])) + 35 * c[3] *
                (24832 * c[4]^2 + 162675 * (3 * c[5]^2)) + 4 * (2097900 *
                c[4]^2 * (2 * c[5]) + 7546875 * (4 * c[5]^3))) + 27 * c[4] *
                (1028300 * c[3]^3 + 50 * c[3]^2 * (594055 * (2 * c[5])) +
                700 * c[3] * (27904 * c[4]^2 + 598575 * (3 * c[5]^2)) + 3 *
              (68908000 * c[4]^2 * (2 * c[5]) + 806378125 * (4 * c[5]^3))))
  j <- 120 * (32 * c[2]^6 + 3456 * c[2]^5 * c[4] + 6 * c[1]^5 * c[5] + 3 *
              c[1]^4 * (9 * c[3]^2 + 16 * c[2] * c[4] + 168 * c[4]^2 + 330 *
              c[3] * c[5] + 2850 * c[5]^2) + 720 * c[2]^4 * (15 * c[3]^2 +
              224 * c[4]^2 + 490 * c[3] * c[5] + 4410 * c[5]^2) + 6048 *
              c[2]^3 * c[4] * (125 * c[3]^2 + 704 * c[4]^2 + 4480 * c[3] *
              c[5] + 44100 * c[5]^2) + 12 * c[1]^3 * (4 * c[2]^2 * (3 * c[3] +
              50 * c[5]) + 60 * c[2] * c[4] * (7 * c[3] + 114 * c[5]) + 3 *
              (24 * c[3]^3 + 1192 * c[3] * c[4]^2 + 1170 * c[3]^2 * c[5] +
                 20440 * c[4]^2 * c[5] +
          20150 * c[3] * c[5]^2 + 124875 * c[5]^3)) + 216 * c[2]^2 * (945 *
              c[3]^4 + 67620 * c[3]^3 * c[5] + 560 * c[3]^2 * (178 * c[4]^2 +
          3555 * c[5]^2) + 315 * c[3] * (12416 * c[4]^2 * c[5] + 90375 *
              c[5]^3) + 6 * (52224 * c[4]^4 + 6993000 * c[4]^2 * c[5]^2 +
          27671875 * c[5]^4)) + 6 * c[1]^2 * (8 * c[2]^4 + 480 * c[2]^3 *
              c[4] + 180 * c[2]^2 * (4 * c[3]^2 + 64 * c[4]^2 + 125 * c[3] *
              c[5] + 1050 * c[5]^2) + 72 * c[2] * c[4] * (327 * c[3]^2 +
              1848 * c[4]^2 + 10955 * c[3] * c[5] + 99900 * c[5]^2) + 9 *
                (225 * c[3]^4 + 22824 * c[3]^2 * c[4]^2 + 69632 * c[4]^4 +
                   15090 * c[3]^3 * c[5] + 830240 * c[3] * c[4]^2 * c[5] +
                   412925 * c[3]^2 * c[5]^2 + 8239800 * c[4]^2 * c[5]^2 +
                   5475750 * c[3] * c[5]^3 + 29636250 * c[5]^4)) + 1296 *
            c[2] * c[4] * (5910 * c[3]^4 + 462735 * c[3]^3 * c[5] + c[3]^2 *
                             (228240 * c[4]^2 + 14851375 * c[5]^2) +
          175 * c[3] * (55808 * c[4]^2 * c[5] + 1316865 * c[5]^3) + 3 *
            (158464 * c[4]^4 + 37899400 * c[4]^2 * c[5]^2 +
          483826875 * c[5]^4)) + 27 * (9945 * c[3]^6 + 92930048 * c[4]^6 +
                                         1166130 * c[3]^5 * c[5] +
          35724729600 * c[4]^4 * c[5]^2 + 977816385000 * c[4]^2 * c[5]^4 +
            1907724656250 * c[5]^6 +
          180 * c[3]^4 * (16082 * c[4]^2 + 345905 * c[5]^2) + 140 * c[3]^3 *
            (1765608 * c[4]^2 * c[5] + 13775375 * c[5]^3) +
          15 * c[3]^2 * (4076032 * c[4]^4 + 574146160 * c[4]^2 * c[5]^2 +
                           2424667875 * c[5]^4) +
          210 * c[3] * (13526272 * c[4]^4 * c[5] + 687499200 * c[4]^2 *
                          c[5]^3 + 1876468125 * c[5]^5)) +
          18 * c[1] * (80 * c[2]^4 * (c[3] + 15 * c[5]) + 160 * c[2]^3 *
                         c[4] * (32 * c[3] + 525 * c[5]) + 12 * c[2]^2 *
                         (225 * c[3]^3 +
          11088 * c[3] * c[4]^2 + 11130 * c[3]^2 * c[5] + 199360 * c[4]^2 *
            c[5] + 200900 * c[3] * c[5]^2 + 1323000 * c[5]^3) +
          24 * c[2] * c[4] * (3885 * c[3]^3 + 69632 * c[3] * c[4]^2 +
                                209800 * c[3]^2 * c[5] + 1369440 * c[4]^2 *
          c[5] + 4134375 * c[3] * c[5]^2 + 29636250 * c[5]^3) + 9 * (540 *
          c[3]^5 + 48585 * c[3]^4 * c[5] + 20 * c[3]^3 * (4856 * c[4]^2 +
          95655 * c[5]^2) + 80 * c[3]^2 * (71597 * c[4]^2 * c[5] + 513625 *
                                             c[5]^3) +
          4 * c[3] * (237696 * c[4]^4 + 30726500 * c[4]^2 * c[5]^2 +
                        119844375 * c[5]^4) +
          5 * c[5] * (4076032 * c[4]^4 + 191074800 * c[4]^2 * c[5]^2 +
                        483826875 * c[5]^4))))
  j_c1 <- 120 * (6 * (5 * c[1]^4) * c[5] + 3 * (4 * c[1]^3) * (9 * c[3]^2 +
                  16 * c[2] * c[4] + 168 * c[4]^2 + 330 * c[3] * c[5] +
              2850 * c[5]^2) + 12 * (3 * c[1]^2) * (4 * c[2]^2 * (3 * c[3] +
                                                    50 * c[5]) +
              60 * c[2] * c[4] * (7 * c[3] + 114 * c[5]) + 3 * (24 * c[3]^3 +
                1192 * c[3] * c[4]^2 + 1170 * c[3]^2 * c[5] +
              20440 * c[4]^2 * c[5] + 20150 * c[3] * c[5]^2 + 124875 *
                c[5]^3)) + 6 * (2 * c[1]) * (8 * c[2]^4 +
              480 * c[2]^3 * c[4] + 180 * c[2]^2 * (4 * c[3]^2 + 64 * c[4]^2 +
                                125 * c[3] * c[5] + 1050 * c[5]^2) +
                72 * c[2] * c[4] * (327 * c[3]^2 + 1848 * c[4]^2 + 10955 *
                c[3] * c[5] + 99900 * c[5]^2) + 9 * (225 * c[3]^4 +
              22824 * c[3]^2 * c[4]^2 + 69632 * c[4]^4 + 15090 * c[3]^3 *
                c[5] + 830240 * c[3] * c[4]^2 * c[5] +
                412925 * c[3]^2 * c[5]^2 + 8239800 * c[4]^2 * c[5]^2 +
                5475750 * c[3] * c[5]^3 + 29636250 * c[5]^4)) +
              18 * (80 * c[2]^4 * (c[3] + 15 * c[5]) + 160 * c[2]^3 * c[4] *
                      (32 * c[3] + 525 * c[5]) +
              12 * c[2]^2 * (225 * c[3]^3 + 11088 * c[3] * c[4]^2 + 11130 *
                               c[3]^2 * c[5] + 199360 * c[4]^2 * c[5] +
              200900 * c[3] * c[5]^2 + 1323000 * c[5]^3) + 24 * c[2] * c[4] *
                (3885 * c[3]^3 + 69632 * c[3] * c[4]^2 +
              209800 * c[3]^2 * c[5] + 1369440 * c[4]^2 * c[5] + 4134375 *
                c[3] * c[5]^2 + 29636250 * c[5]^3) + 9 * (540 * c[3]^5 +
              48585 * c[3]^4 * c[5] + 20 * c[3]^3 * (4856 * c[4]^2 + 95655 *
                c[5]^2) + 80 * c[3]^2 * (71597 * c[4]^2 * c[5] +
              513625 * c[5]^3) + 4 * c[3] * (237696 * c[4]^4 + 30726500 *
                c[4]^2 * c[5]^2 + 119844375 * c[5]^4) +
              5 * c[5] * (4076032 * c[4]^4 + 191074800 * c[4]^2 * c[5]^2 +
                            483826875 * c[5]^4))))
  j_c2 <- 120 * (32 * (6 * c[2]^5) + 3456 * (5 * c[2]^4) * c[4] + 3 * c[1]^4 *
                   (16 * c[4]) + 720 * (4 * c[2]^3) * (15 * c[3]^2 +
              224 * c[4]^2 + 490 * c[3] * c[5] + 4410 * c[5]^2) + 6048 *
                (3 * c[2]^2) * c[4] * (125 * c[3]^2 +
              704 * c[4]^2 + 4480 * c[3] * c[5] + 44100 * c[5]^2) + 12 *
                c[1]^3 * (4 * (2 * c[2]) * (3 * c[3] + 50 * c[5]) +
              60 * c[4] * (7 * c[3] + 114 * c[5])) + 216 * (2 * c[2]) *
                (945 * c[3]^4 + 67620 * c[3]^3 * c[5] +
              560 * c[3]^2 * (178 * c[4]^2 + 3555 * c[5]^2) + 315 * c[3] *
                (12416 * c[4]^2 * c[5] + 90375 * c[5]^3) +
                6 * (52224 * c[4]^4 + 6993000 * c[4]^2 * c[5]^2 + 27671875 *
                       c[5]^4)) +
              6 * c[1]^2 * (8 * (4 * c[2]^3) + 480 * (3 * c[2]^2) * c[4] +
                              180 * (2 * c[2]) * (4 * c[3]^2 +
              64 * c[4]^2 + 125 * c[3] * c[5] + 1050 * c[5]^2) + 72 * c[4] *
                (327 * c[3]^2 + 1848 * c[4]^2 + 10955 * c[3] * c[5] +
              99900 * c[5]^2)) + 1296 * c[4] * (5910 * c[3]^4 + 462735 *
              c[3]^3 * c[5] + c[3]^2 * (228240 * c[4]^2 + 14851375 * c[5]^2) +
              175 * c[3] * (55808 * c[4]^2 * c[5] + 1316865 * c[5]^3) + 3 *
                (158464 * c[4]^4 + 37899400 * c[4]^2 * c[5]^2 +
              483826875 * c[5]^4)) + 18 * c[1] * (80 * (4 * c[2]^3) * (c[3] +
              15 * c[5]) + 160 * (3 * c[2]^2) * c[4] * (32 * c[3] + 525 *
                                                          c[5]) +
              12 * (2 * c[2]) * (225 * c[3]^3 + 11088 * c[3] * c[4]^2 +
              11130 * c[3]^2 * c[5] + 199360 * c[4]^2 * c[5] + 200900 * c[3] *
                c[5]^2 +
              1323000 * c[5]^3) + 24 * c[4] * (3885 * c[3]^3 + 69632 * c[3] *
                                                 c[4]^2 +
              209800 * c[3]^2 * c[5] + 1369440 * c[4]^2 * c[5] + 4134375 *
                c[3] * c[5]^2 + 29636250 * c[5]^3)))
  j_c3 <- 120 * (3 * c[1]^4 * (9 * (2 * c[3]) + 330 * c[5]) + 720 * c[2]^4 *
                   (15 * (2 * c[3]) + 490 * c[5]) +
                6048 * c[2]^3 * c[4] * (125 * (2 * c[3]) + 4480 * c[5]) +
                  12 * c[1]^3 * (4 * c[2]^2 * 3 + 60 * c[2] * c[4] * 7 +
              3 * (24 * (3 * c[3]^2) + 1192 * c[4]^2 + 1170 * (2 * c[3]) *
              c[5] + 20150 * c[5]^2)) + 216 * c[2]^2 * (945 * (4 * c[3]^3) +
            67620 * (3 * c[3]^2) * c[5] + 560 * (2 * c[3]) * (178 * c[4]^2 +
                                                            3555 * c[5]^2) +
            315 * (12416 * c[4]^2 * c[5] + 90375 * c[5]^3)) + 6 * c[1]^2 *
              (180 * c[2]^2 * (4 * (2 * c[3]) + 125 * c[5]) +
          72 * c[2] * c[4] * (327 * (2 * c[3]) + 10955 * c[5]) + 9 * (225 *
              (4 * c[3]^3) + 22824 * (2 * c[3]) * c[4]^2 +
          15090 * (3 * c[3]^2) * c[5] + 830240 * c[4]^2 * c[5] + 412925 *
            (2 * c[3]) * c[5]^2 + 5475750 * c[5]^3)) +
            1296 * c[2] * c[4] * (5910 * (4 * c[3]^3) + 462735 *
            (3 * c[3]^2) * c[5] + 2 * c[3] * (228240 * c[4]^2 + 14851375 *
                                                c[5]^2) +
          175 * (55808 * c[4]^2 * c[5] + 1316865 * c[5]^3)) + 27 * (9945 *
          (6 * c[3]^5) + 1166130 * (5 * c[3]^4) * c[5] + 180 * (4 * c[3]^3) *
            (16082 * c[4]^2 +
          345905 * c[5]^2) + 140 * (3 * c[3]^2) * (1765608 * c[4]^2 * c[5] +
            13775375 * c[5]^3) + 15 * (2 * c[3]) * (4076032 * c[4]^4 +
          574146160 * c[4]^2 * c[5]^2 + 2424667875 * c[5]^4) + 210 *
            (13526272 * c[4]^4 * c[5] + 687499200 * c[4]^2 * c[5]^3 +
          1876468125 * c[5]^5)) + 18 * c[1] * (80 * c[2]^4 + 160 * c[2]^3 *
            c[4] * 32 + 12 * c[2]^2 * (225 * (3 * c[3]^2) + 11088 * c[4]^2 +
          11130 * (2 * c[3]) * c[5] + 200900 * c[5]^2) + 24 * c[2] * c[4] *
          (3885 * (3 * c[3]^2) + 69632 * c[4]^2 + 209800 * (2 * c[3]) * c[5] +
          4134375 * c[5]^2) + 9 * (540 * (5 * c[3]^4) + 48585 * (4 * c[3]^3) *
                c[5] + 20 * (3 * c[3]^2) * (4856 * c[4]^2 + 95655 * c[5]^2) +
          80 * (2 * c[3]) * (71597 * c[4]^2 * c[5] + 513625 * c[5]^3) + 4 *
          (237696 * c[4]^4 + 30726500 * c[4]^2 * c[5]^2 + 119844375 *
             c[5]^4))))
  j_c4 <- 120 * (3456 * c[2]^5 + 3 * c[1]^4 * (16 * c[2] + 168 * (2 * c[4])) +
                   720 * c[2]^4 * (224 * (2 * c[4])) + (6048 * c[2]^3 * (125 *
                  c[3]^2 + 704 * c[4]^2 + 4480 * c[3] * c[5] + 44100 *
                    c[5]^2) + 6048 * c[2]^3 * c[4] * (704 * (2 * c[4]))) +
              12 * c[1]^3 * (60 * c[2] * (7 * c[3] + 114 * c[5]) + 3 *
              (1192 * c[3] * (2 * c[4]) + 20440 * (2 * c[4]) * c[5])) + 216 *
                c[2]^2 * (560 * c[3]^2 * (178 * (2 * c[4])) +
              315 * c[3] * (12416 * (2 * c[4]) * c[5]) + 6 * (52224 * (4 *
                c[4]^3) + 6993000 * (2 * c[4]) * c[5]^2)) + 6 * c[1]^2 *
                (480 * c[2]^3 + 180 * c[2]^2 * (64 * (2 * c[4])) + (72 *
                c[2] * (327 * c[3]^2 + 1848 * c[4]^2 + 10955 * c[3] * c[5] +
                          99900 * c[5]^2) +
              72 * c[2] * c[4] * (1848 * (2 * c[4]))) + 9 * (22824 * c[3]^2 *
                (2 * c[4]) + 69632 * (4 * c[4]^3) + 830240 * c[3] *
                  (2 * c[4]) * c[5] +
              8239800 * (2 * c[4]) * c[5]^2)) + (1296 * c[2] *
                (5910 * c[3]^4 + 462735 * c[3]^3 * c[5] + c[3]^2 *
                   (228240 * c[4]^2 +
              14851375 * c[5]^2) + 175 * c[3] * (55808 * c[4]^2 * c[5] +
                1316865 * c[5]^3) + 3 * (158464 * c[4]^4 +
              37899400 * c[4]^2 * c[5]^2 + 483826875 * c[5]^4)) + 1296 *
                c[2] * c[4] * (c[3]^2 * (228240 * (2 * c[4])) +
              175 * c[3] * (55808 * (2 * c[4]) * c[5]) + 3 * (158464 *
                (4 * c[4]^3) + 37899400 * (2 * c[4]) * c[5]^2))) +
                27 * (92930048 * (6 * c[4]^5) + 35724729600 * (4 * c[4]^3) *
                c[5]^2 + 977816385000 * (2 * c[4]) * c[5]^4 + 180 * c[3]^4 *
                  (16082 * (2 * c[4])) +
                140 * c[3]^3 * (1765608 * (2 * c[4]) * c[5]) + 15 * c[3]^2 *
                  (4076032 * (4 * c[4]^3) + 574146160 * (2 * c[4]) * c[5]^2) +
                210 * c[3] * (13526272 * (4 * c[4]^3) * c[5] + 687499200 *
                                (2 * c[4]) * c[5]^3)) +
              18 * c[1] * (160 * c[2]^3 * (32 * c[3] + 525 * c[5]) + 12 *
                             c[2]^2 * (11088 * c[3] * (2 * c[4]) +
              199360 * (2 * c[4]) * c[5]) + (24 * c[2] * (3885 * c[3]^3 +
                  69632 * c[3] * c[4]^2 + 209800 * c[3]^2 * c[5] +
              1369440 * c[4]^2 * c[5] + 4134375 * c[3] * c[5]^2 + 29636250 *
                c[5]^3) + 24 * c[2] * c[4] * (69632 * c[3] * (2 * c[4]) +
              1369440 * (2 * c[4]) * c[5])) + 9 * (20 * c[3]^3 * (4856 *
                (2 * c[4])) + 80 * c[3]^2 * (71597 * (2 * c[4]) * c[5]) +
            4 * c[3] * (237696 * (4 * c[4]^3) + 30726500 * (2 * c[4]) *
                c[5]^2) + 5 * c[5] * (4076032 * (4 * c[4]^3) + 191074800 *
                                        (2 * c[4]) * c[5]^2))))
  j_c5 <- 120 * (6 * c[1]^5 + 3 * c[1]^4 * (330 * c[3] + 2850 * (2 * c[5])) +
                   720 * c[2]^4 * (490 * c[3] + 4410 * (2 * c[5])) +
            6048 * c[2]^3 * c[4] * (4480 * c[3] + 44100 * (2 * c[5])) +
              12 * c[1]^3 * (4 * c[2]^2 * 50 +
            60 * c[2] * c[4] * 114 + 3 * (1170 * c[3]^2 + 20440 * c[4]^2 +
              20150 * c[3] * (2 * c[5]) + 124875 * (3 * c[5]^2))) +
            216 * c[2]^2 * (67620 * c[3]^3 + 560 * c[3]^2 *
                              (3555 * (2 * c[5])) +
            315 * c[3] * (12416 * c[4]^2 + 90375 * (3 * c[5]^2)) + 6 *
              (6993000 * c[4]^2 * (2 * c[5]) + 27671875 * (4 * c[5]^3))) +
            6 * c[1]^2 * (180 * c[2]^2 * (125 * c[3] + 1050 * (2 * c[5])) +
              72 * c[2] * c[4] * (10955 * c[3] + 99900 * (2 * c[5])) +
            9 * (15090 * c[3]^3 + 830240 * c[3] * c[4]^2 + 412925 * c[3]^2 *
              (2 * c[5]) + 8239800 * c[4]^2 * (2 * c[5]) + 5475750 * c[3] *
                (3 * c[5]^2) +
            29636250 * (4 * c[5]^3))) + 1296 * c[2] * c[4] *
              (462735 * c[3]^3 + c[3]^2 * (14851375 * (2 * c[5])) +
                 175 * c[3] * (55808 * c[4]^2 +
            1316865 * (3 * c[5]^2)) + 3 * (37899400 * c[4]^2 * (2 * c[5]) +
              483826875 * (4 * c[5]^3))) + 27 * (1166130 * c[3]^5 +
            35724729600 * c[4]^4 * (2 * c[5]) + 977816385000 * c[4]^2 *
              (4 * c[5]^3) + 1907724656250 * (6 * c[5]^5) +
            180 * c[3]^4 * (345905 * (2 * c[5])) + 140 * c[3]^3 * (1765608 *
              c[4]^2 + 13775375 * (3 * c[5]^2)) +
            15 * c[3]^2 * (574146160 * c[4]^2 * (2 * c[5]) + 2424667875 *
              (4 * c[5]^3)) + 210 * c[3] * (13526272 * c[4]^4 +
            687499200 * c[4]^2 * (3 * c[5]^2) + 1876468125 * (5 * c[5]^4))) +
              18 * c[1] * (80 * c[2]^4 * 15 +
            160 * c[2]^3 * c[4] * 525 + 12 * c[2]^2 * (11130 * c[3]^2 +
              199360 * c[4]^2 + 200900 * c[3] * (2 * c[5]) +
            1323000 * (3 * c[5]^2)) + 24 * c[2] * c[4] * (209800 * c[3]^2 +
              1369440 * c[4]^2 + 4134375 * c[3] * (2 * c[5]) +
                29636250 * (3 * c[5]^2)) +
            9 * (48585 * c[3]^4 + 20 * c[3]^3 * (95655 * (2 * c[5])) +
                   80 * c[3]^2 * (71597 * c[4]^2 + 513625 * (3 * c[5]^2)) +
            4 * c[3] * (30726500 * c[4]^2 * (2 * c[5]) + 119844375 *
                          (4 * c[5]^3)) + (5 * (4076032 * c[4]^4 +
            191074800 * c[4]^2 * c[5]^2 + 483826875 * c[5]^4) + 5 * c[5] *
              (191074800 * c[4]^2 * (2 * c[5]) + 483826875 * (4 * c[5]^3))))))
  F[1] <- 1 - g
  F[2] <- g1 - h
  F[3] <- g3 - i
  F[4] <- g4 - j
  F[5] <- f_c1 - c[6] * g_c1 - c[7] * h_c1 - c[8] * i_c1 - c[9] * j_c1
  F[6] <- f_c2 - c[6] * g_c2 - c[7] * h_c2 - c[8] * i_c2 - c[9] * j_c2
  F[7] <- f_c3 - c[6] * g_c3 - c[7] * h_c3 - c[8] * i_c3 - c[9] * j_c3
  F[8] <- f_c4 - c[6] * g_c4 - c[7] * h_c4 - c[8] * i_c4 - c[9] * j_c4
  F[9] <- f_c5 - c[6] * g_c5 - c[7] * h_c5 - c[8] * i_c5 - c[9] * j_c5
  return(F)
}