R Under development (unstable) (2023-12-04 r85659 ucrt) -- "Unsuffered Consequences"
Copyright (C) 2023 The R Foundation for Statistical Computing
Platform: x86_64-w64-mingw32/x64

R is free software and comes with ABSOLUTELY NO WARRANTY.
You are welcome to redistribute it under certain conditions.
Type 'license()' or 'licence()' for distribution details.

R is a collaborative project with many contributors.
Type 'contributors()' for more information and
'citation()' on how to cite R or R packages in publications.

Type 'demo()' for some demos, 'help()' for on-line help, or
'help.start()' for an HTML browser interface to help.
Type 'q()' to quit R.

> ## Copyright (C) 2012 Marius Hofert, Ivan Kojadinovic, Martin Maechler, and Jun Yan
> ##
> ## This program is free software; you can redistribute it and/or modify it under
> ## the terms of the GNU General Public License as published by the Free Software
> ## Foundation; either version 3 of the License, or (at your option) any later
> ## version.
> ##
> ## This program is distributed in the hope that it will be useful, but WITHOUT
> ## ANY WARRANTY; without even the implied warranty of MERCHANTABILITY or FITNESS
> ## FOR A PARTICULAR PURPOSE. See the GNU General Public License for more
> ## details.
> ##
> ## You should have received a copy of the GNU General Public License along with
> ## this program; if not, see <http://www.gnu.org/licenses/>.
> 
> ### pdf = dCopula()  ---  but then also pCopula()
> ### ===   ~~~~~~~~~                     ~~~~~~~~~
> 
> require(copula)
Loading required package: copula
> ## library(fCopulae)
> 
> source(system.file("Rsource", "utils.R", package="copula", mustWork=TRUE))
Loading required package: tools
> ## tryCatch.W.E() etc
> ## All non-virtual copula classes: ../inst/Rsource/cops.R
> source(system.file("Rsource", "cops.R", package="copula", mustWork=TRUE))
> ## --> copcl, copObs, copBnds,  excl.2 , copO.2, copBnd.2
> (doExtras <- copula:::doExtras())
[1] FALSE
> 
> 
> ### preparation for a grid
> n1 <- 17
> n2 <- 21
> eps <- .001 ## <- not going to very extremes
> u1 <- seq(0, 1, length=n1); u1[1] <- eps; u1[n1] <- 1-eps
> u2 <- seq(0, 1, length=n2); u2[1] <- eps; u2[n2] <- 1-eps
> 
> ### d=2 ########################################################################
> 
> Exp.grid <- function(...)
+     as.matrix(expand.grid(..., KEEP.OUT.ATTRS = FALSE))
> 
> umat <- Exp.grid(u1=u1, u2=u2)
> 
> ## all copulas give the same tau except amh (TODO: *range* of tau's ??)
> tau <- 0.5
> 
> if(!dev.interactive(orNone=TRUE)) pdf("densCop_2d.pdf")
> 
> fCols <- colorRampPalette(c("red", "white", "blue"), space = "Lab")
> 
> options(width = 137)# -> nicer table printing
> showProc.time()
Time (user system elapsed): 0.09 0 0.1 
> 
> 
> ## frankCopula
> theta.fr <- iTau(frankCopula(), tau)
> dcop <- matrix(dCopula(umat, frankCopula(param=theta.fr, dim = 2)),
+                n1, n2)
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( frankCopula(%.4g) )", theta.fr))
> round(dcop, 3)
       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 5.690 4.307 3.240 2.436 1.831 1.376 1.033 0.776 0.583 0.438 0.329 0.247 0.185 0.139 0.104 0.078 0.059 0.044 0.033 0.025 0.019
 [2,] 4.012 3.530 3.007 2.495 2.027 1.617 1.273 0.991 0.764 0.586 0.447 0.340 0.257 0.194 0.147 0.111 0.083 0.063 0.047 0.035 0.027
 [3,] 2.810 2.774 2.630 2.400 2.114 1.806 1.502 1.223 0.977 0.770 0.600 0.463 0.355 0.271 0.206 0.156 0.118 0.089 0.067 0.051 0.038
 [4,] 1.967 2.113 2.185 2.170 2.069 1.898 1.680 1.441 1.203 0.982 0.787 0.621 0.484 0.375 0.288 0.220 0.167 0.127 0.096 0.072 0.055
 [5,] 1.376 1.571 1.740 1.856 1.902 1.871 1.768 1.608 1.412 1.203 0.999 0.812 0.648 0.510 0.397 0.306 0.235 0.179 0.136 0.103 0.078
 [6,] 0.962 1.148 1.339 1.514 1.651 1.732 1.745 1.688 1.569 1.406 1.219 1.027 0.844 0.681 0.540 0.423 0.328 0.252 0.193 0.147 0.112
 [7,] 0.673 0.828 1.004 1.187 1.363 1.513 1.618 1.663 1.640 1.554 1.417 1.247 1.065 0.885 0.720 0.576 0.453 0.353 0.273 0.209 0.160
 [8,] 0.470 0.592 0.738 0.903 1.079 1.256 1.415 1.538 1.608 1.615 1.556 1.442 1.288 1.114 0.936 0.768 0.618 0.490 0.383 0.297 0.230
 [9,] 0.329 0.420 0.535 0.671 0.827 0.999 1.176 1.343 1.482 1.574 1.607 1.574 1.482 1.343 1.176 0.999 0.827 0.671 0.535 0.420 0.329
[10,] 0.230 0.297 0.383 0.490 0.618 0.768 0.936 1.114 1.288 1.442 1.556 1.615 1.608 1.538 1.415 1.256 1.079 0.903 0.738 0.592 0.470
[11,] 0.160 0.209 0.273 0.353 0.453 0.576 0.720 0.885 1.065 1.247 1.417 1.554 1.640 1.663 1.618 1.513 1.363 1.187 1.004 0.828 0.673
[12,] 0.112 0.147 0.193 0.252 0.328 0.423 0.540 0.681 0.844 1.027 1.219 1.406 1.569 1.688 1.745 1.732 1.651 1.514 1.339 1.148 0.962
[13,] 0.078 0.103 0.136 0.179 0.235 0.306 0.397 0.510 0.648 0.812 0.999 1.203 1.412 1.608 1.768 1.871 1.902 1.856 1.740 1.571 1.376
[14,] 0.055 0.072 0.096 0.127 0.167 0.220 0.288 0.375 0.484 0.621 0.787 0.982 1.203 1.441 1.680 1.898 2.069 2.170 2.185 2.113 1.967
[15,] 0.038 0.051 0.067 0.089 0.118 0.156 0.206 0.271 0.355 0.463 0.600 0.770 0.977 1.223 1.502 1.806 2.114 2.400 2.630 2.774 2.810
[16,] 0.027 0.035 0.047 0.063 0.083 0.111 0.147 0.194 0.257 0.340 0.447 0.586 0.764 0.991 1.273 1.617 2.027 2.495 3.007 3.530 4.012
[17,] 0.019 0.025 0.033 0.044 0.059 0.078 0.104 0.139 0.185 0.247 0.329 0.438 0.583 0.776 1.033 1.376 1.831 2.436 3.240 4.307 5.690
> 
> 
> ## claytonCopula
> (theta.cl <- iTau(claytonCopula(), tau))
[1] 2
> stopifnot(all.equal(theta.cl, copClayton@iTau(tau), tolerance = 1e-13))
> dcop <- matrix(dCopula(umat, claytonCopula(param=theta.cl, dim = 2)),
+                n1, n2)
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf(" Density( claytonCopula(%.4g) )", theta.cl))
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( claytonCopula(%.4g) )", theta.cl))
> round(dcop, 3)
         [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 530.331 0.024 0.003 0.001 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.000
 [2,]   0.012 8.953 5.175 2.346 1.171 0.650 0.394 0.255 0.174 0.124 0.091 0.069 0.053 0.042 0.034 0.028 0.023 0.019 0.016 0.014 0.012
 [3,]   0.002 2.664 4.528 3.803 2.643 1.772 1.203 0.839 0.601 0.443 0.334 0.258 0.203 0.162 0.131 0.108 0.090 0.075 0.064 0.054 0.047
 [4,]   0.000 0.964 2.482 3.077 2.856 2.341 1.826 1.403 1.079 0.837 0.657 0.522 0.420 0.342 0.281 0.234 0.196 0.166 0.142 0.122 0.106
 [5,]   0.000 0.438 1.354 2.088 2.372 2.297 2.041 1.734 1.441 1.187 0.976 0.805 0.667 0.556 0.467 0.395 0.336 0.288 0.248 0.215 0.188
 [6,]   0.000 0.232 0.788 1.379 1.791 1.966 1.949 1.815 1.627 1.425 1.233 1.060 0.909 0.780 0.671 0.578 0.501 0.435 0.380 0.333 0.294
 [7,]   0.000 0.137 0.490 0.928 1.317 1.584 1.712 1.724 1.655 1.538 1.400 1.257 1.118 0.990 0.874 0.771 0.681 0.602 0.534 0.474 0.423
 [8,]   0.000 0.087 0.323 0.642 0.970 1.246 1.441 1.547 1.576 1.547 1.477 1.384 1.278 1.170 1.063 0.962 0.869 0.783 0.706 0.636 0.575
 [9,]   0.000 0.059 0.223 0.459 0.723 0.976 1.188 1.344 1.441 1.483 1.481 1.445 1.385 1.310 1.226 1.140 1.054 0.971 0.892 0.818 0.751
[10,]   0.000 0.042 0.160 0.337 0.548 0.768 0.973 1.148 1.284 1.378 1.432 1.450 1.440 1.407 1.358 1.298 1.231 1.160 1.089 1.018 0.951
[11,]   0.000 0.030 0.118 0.254 0.422 0.609 0.796 0.972 1.126 1.253 1.349 1.414 1.452 1.464 1.456 1.431 1.393 1.345 1.291 1.233 1.173
[12,]   0.000 0.023 0.090 0.195 0.331 0.488 0.654 0.821 0.980 1.123 1.248 1.350 1.429 1.485 1.521 1.537 1.536 1.522 1.495 1.460 1.419
[13,]   0.000 0.018 0.070 0.153 0.263 0.395 0.540 0.694 0.849 0.999 1.140 1.268 1.381 1.477 1.554 1.615 1.658 1.685 1.698 1.698 1.688
[14,]   0.000 0.014 0.055 0.122 0.213 0.323 0.449 0.588 0.734 0.884 1.033 1.178 1.316 1.445 1.562 1.666 1.756 1.833 1.896 1.944 1.980
[15,]   0.000 0.011 0.044 0.099 0.174 0.267 0.377 0.500 0.636 0.780 0.931 1.086 1.241 1.396 1.547 1.693 1.832 1.963 2.085 2.196 2.295
[16,]   0.000 0.009 0.036 0.081 0.144 0.223 0.318 0.428 0.552 0.688 0.836 0.994 1.161 1.335 1.515 1.699 1.886 2.075 2.264 2.451 2.633
[17,]   0.000 0.008 0.030 0.068 0.120 0.188 0.271 0.368 0.481 0.609 0.751 0.909 1.081 1.269 1.471 1.688 1.920 2.166 2.427 2.703 2.988
> 
> 
> 
> ## gumbelCopula
> theta.gu <- iTau(gumbelCopula(), tau)
> stopifnot(all.equal(theta.gu, copGumbel@iTau(tau), tolerance = 1e-13))
> dcop <- matrix(dCopula(umat, gumbelCopula(param=theta.gu, dim = 2)),
+                n1, n2)
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( gumbelCopula(%.4g) )", theta.gu))
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( gumbelCopula(%.4g) )", theta.gu))
> round(dcop, 3)
        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]   [,21]
 [1,] 31.525 4.442 2.348 1.502 1.048 0.768 0.581 0.448 0.351 0.277 0.220 0.174 0.138 0.108 0.083 0.063 0.046 0.032 0.019 0.009   0.000
 [2,]  3.672 3.352 2.734 2.242 1.847 1.526 1.262 1.044 0.863 0.711 0.583 0.475 0.384 0.306 0.240 0.183 0.135 0.093 0.057 0.026   0.000
 [3,]  1.851 2.490 2.364 2.156 1.927 1.699 1.481 1.280 1.095 0.929 0.781 0.649 0.532 0.430 0.340 0.262 0.194 0.134 0.083 0.038   0.001
 [4,]  1.140 1.896 1.989 1.960 1.869 1.741 1.592 1.432 1.269 1.109 0.956 0.811 0.677 0.555 0.444 0.345 0.257 0.179 0.111 0.051   0.001
 [5,]  0.768 1.464 1.650 1.728 1.739 1.701 1.626 1.522 1.399 1.262 1.117 0.971 0.828 0.690 0.560 0.439 0.330 0.231 0.144 0.067   0.001
 [6,]  0.543 1.139 1.355 1.489 1.568 1.602 1.596 1.555 1.482 1.383 1.263 1.129 0.985 0.838 0.692 0.550 0.418 0.295 0.184 0.086   0.002
 [7,]  0.396 0.888 1.101 1.258 1.377 1.462 1.514 1.532 1.516 1.467 1.386 1.278 1.147 1.000 0.844 0.684 0.526 0.376 0.236 0.111   0.002
 [8,]  0.294 0.693 0.887 1.045 1.181 1.297 1.390 1.459 1.498 1.504 1.475 1.409 1.308 1.176 1.019 0.845 0.662 0.479 0.305 0.143   0.003
 [9,]  0.220 0.538 0.706 0.853 0.989 1.117 1.236 1.340 1.426 1.487 1.516 1.507 1.454 1.357 1.217 1.040 0.836 0.617 0.397 0.189   0.004
[10,]  0.164 0.414 0.554 0.682 0.808 0.935 1.062 1.187 1.305 1.412 1.498 1.553 1.567 1.529 1.432 1.273 1.059 0.804 0.528 0.254   0.005
[11,]  0.122 0.314 0.426 0.533 0.643 0.758 0.880 1.009 1.143 1.280 1.411 1.530 1.620 1.665 1.644 1.541 1.346 1.065 0.721 0.353   0.007
[12,]  0.089 0.233 0.320 0.404 0.494 0.591 0.699 0.819 0.952 1.098 1.256 1.421 1.582 1.723 1.814 1.821 1.702 1.432 1.019 0.514   0.010
[13,]  0.063 0.167 0.231 0.294 0.363 0.439 0.527 0.627 0.744 0.881 1.040 1.223 1.430 1.654 1.873 2.046 2.101 1.944 1.502 0.801   0.016
[14,]  0.042 0.112 0.156 0.201 0.249 0.304 0.368 0.444 0.535 0.646 0.782 0.950 1.160 1.419 1.732 2.086 2.421 2.587 2.315 1.382   0.028
[15,]  0.025 0.068 0.094 0.122 0.152 0.186 0.227 0.276 0.336 0.411 0.507 0.631 0.795 1.019 1.328 1.760 2.352 3.075 3.583 2.789   0.064
[16,]  0.011 0.031 0.043 0.055 0.069 0.085 0.104 0.127 0.156 0.192 0.239 0.302 0.389 0.513 0.698 0.993 1.497 2.428 4.243 6.613   0.256
[17,]  0.000 0.000 0.001 0.001 0.001 0.001 0.002 0.002 0.002 0.003 0.004 0.004 0.006 0.008 0.011 0.016 0.025 0.044 0.100 0.400 354.084
> 
> 
> ## normalCopula
> uB <- cbind(c(0:1, .5), (1:9)/10) # values at boundaries
> fC <- dCopula(uB, normalCopula(0.55))
> stopifnot(is.finite(fC), length(fC)==nrow(uB), fC[-3*(1:3)] == 0)
> theta.n <- iTau(normalCopula(), tau)
> dcop <- matrix(dCopula(umat, normalCopula(param=theta.n, dim = 2)),
+                n1, n2)
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( normalCopula(%.4g) )", theta.n))
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( normalCopula(%.4g) )", theta.n))
> round(dcop, 3)
        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]  [,21]
 [1,] 73.856 4.086 1.421 0.647 0.331 0.181 0.103 0.060 0.035 0.021 0.012 0.007 0.004 0.002 0.001 0.000 0.000 0.000 0.000 0.000  0.000
 [2,]  3.003 3.998 3.093 2.414 1.900 1.500 1.185 0.934 0.732 0.568 0.436 0.329 0.244 0.175 0.122 0.080 0.049 0.027 0.012 0.003  0.000
 [3,]  0.939 2.740 2.582 2.302 2.014 1.741 1.493 1.268 1.067 0.888 0.730 0.590 0.468 0.362 0.271 0.194 0.130 0.079 0.040 0.013  0.000
 [4,]  0.389 1.942 2.095 2.047 1.925 1.771 1.606 1.437 1.270 1.108 0.954 0.809 0.672 0.546 0.431 0.326 0.233 0.152 0.084 0.031  0.000
 [5,]  0.181 1.398 1.683 1.769 1.764 1.707 1.619 1.510 1.389 1.260 1.126 0.991 0.857 0.724 0.595 0.472 0.354 0.245 0.146 0.061  0.000
 [6,]  0.090 1.011 1.339 1.501 1.576 1.593 1.572 1.521 1.448 1.358 1.255 1.142 1.020 0.893 0.761 0.627 0.492 0.358 0.228 0.104  0.001
 [7,]  0.046 0.729 1.054 1.253 1.378 1.451 1.484 1.485 1.459 1.411 1.344 1.260 1.161 1.049 0.925 0.790 0.646 0.492 0.332 0.166  0.003
 [8,]  0.023 0.521 0.817 1.028 1.182 1.293 1.368 1.413 1.431 1.425 1.397 1.348 1.279 1.190 1.083 0.958 0.813 0.648 0.462 0.250  0.006
 [9,]  0.012 0.366 0.622 0.827 0.992 1.126 1.233 1.313 1.370 1.403 1.414 1.403 1.370 1.313 1.233 1.126 0.992 0.827 0.622 0.366  0.012
[10,]  0.006 0.250 0.462 0.648 0.813 0.958 1.083 1.190 1.279 1.348 1.397 1.425 1.431 1.413 1.368 1.293 1.182 1.028 0.817 0.521  0.023
[11,]  0.003 0.166 0.332 0.492 0.646 0.790 0.925 1.049 1.161 1.260 1.344 1.411 1.459 1.485 1.484 1.451 1.378 1.253 1.054 0.729  0.046
[12,]  0.001 0.104 0.228 0.358 0.492 0.627 0.761 0.893 1.020 1.142 1.255 1.358 1.448 1.521 1.572 1.593 1.576 1.501 1.339 1.011  0.090
[13,]  0.000 0.061 0.146 0.245 0.354 0.472 0.595 0.724 0.857 0.991 1.126 1.260 1.389 1.510 1.619 1.707 1.764 1.769 1.683 1.398  0.181
[14,]  0.000 0.031 0.084 0.152 0.233 0.326 0.431 0.546 0.672 0.809 0.954 1.108 1.270 1.437 1.606 1.771 1.925 2.047 2.095 1.942  0.389
[15,]  0.000 0.013 0.040 0.079 0.130 0.194 0.271 0.362 0.468 0.590 0.730 0.888 1.067 1.268 1.493 1.741 2.014 2.302 2.582 2.740  0.939
[16,]  0.000 0.003 0.012 0.027 0.049 0.080 0.122 0.175 0.244 0.329 0.436 0.568 0.732 0.934 1.185 1.500 1.900 2.414 3.093 3.998  3.003
[17,]  0.000 0.000 0.000 0.000 0.000 0.000 0.001 0.002 0.004 0.007 0.012 0.021 0.035 0.060 0.103 0.181 0.331 0.647 1.421 4.086 73.856
> 
> 
> ## tCopula
> fC <- dCopula(uB, tCopula(0.55))
> stopifnot(is.finite(fC), length(fC)==nrow(uB), fC[-3*(1:3)] == 0)
> (theta.t. <- iTau(tCopula(df=10), tau))
[1] 0.7071068
> dcop <- matrix(dCopula(umat, tCopula(param=theta.t., dim = 2, df=10)),
+                n1, n2)
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( tCopula(%.4g, df = 10) )", theta.t.))
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( tCopula(%.4g, df = 10) )", theta.t.))
> round(dcop, 3)
         [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]   [,21]
 [1,] 118.692 2.538 0.869 0.443 0.269 0.180 0.129 0.096 0.074 0.059 0.048 0.039 0.033 0.028 0.023 0.020 0.017 0.015 0.013 0.012   0.021
 [2,]   1.817 4.375 3.283 2.425 1.815 1.377 1.059 0.822 0.644 0.507 0.400 0.317 0.251 0.198 0.156 0.121 0.093 0.070 0.050 0.033   0.012
 [3,]   0.602 2.762 2.750 2.442 2.088 1.754 1.460 1.208 0.995 0.817 0.667 0.542 0.437 0.349 0.276 0.215 0.164 0.121 0.084 0.051   0.014
 [4,]   0.302 1.818 2.145 2.158 2.034 1.852 1.648 1.443 1.247 1.066 0.902 0.755 0.624 0.510 0.409 0.322 0.247 0.182 0.125 0.074   0.017
 [5,]   0.180 1.245 1.648 1.817 1.850 1.800 1.697 1.564 1.415 1.259 1.104 0.954 0.811 0.678 0.556 0.446 0.346 0.256 0.176 0.103   0.020
 [6,]   0.120 0.878 1.262 1.494 1.619 1.666 1.654 1.597 1.507 1.395 1.267 1.131 0.992 0.852 0.716 0.586 0.463 0.348 0.241 0.141   0.024
 [7,]   0.084 0.633 0.966 1.209 1.380 1.491 1.548 1.560 1.533 1.474 1.388 1.281 1.159 1.026 0.886 0.743 0.601 0.460 0.323 0.190   0.030
 [8,]   0.062 0.463 0.739 0.967 1.153 1.298 1.404 1.471 1.502 1.498 1.462 1.397 1.306 1.193 1.062 0.917 0.761 0.597 0.429 0.255   0.037
 [9,]   0.048 0.343 0.564 0.765 0.945 1.104 1.238 1.345 1.423 1.471 1.487 1.471 1.423 1.345 1.238 1.104 0.945 0.765 0.564 0.343   0.048
[10,]   0.037 0.255 0.429 0.597 0.761 0.917 1.062 1.193 1.306 1.397 1.462 1.498 1.502 1.471 1.404 1.298 1.153 0.967 0.739 0.463   0.062
[11,]   0.030 0.190 0.323 0.460 0.601 0.743 0.886 1.026 1.159 1.281 1.388 1.474 1.533 1.560 1.548 1.491 1.380 1.209 0.966 0.633   0.084
[12,]   0.024 0.141 0.241 0.348 0.463 0.586 0.716 0.852 0.992 1.131 1.267 1.395 1.507 1.597 1.654 1.666 1.619 1.494 1.262 0.878   0.120
[13,]   0.020 0.103 0.176 0.256 0.346 0.446 0.556 0.678 0.811 0.954 1.104 1.259 1.415 1.564 1.697 1.800 1.850 1.817 1.648 1.245   0.180
[14,]   0.017 0.074 0.125 0.182 0.247 0.322 0.409 0.510 0.624 0.755 0.902 1.066 1.247 1.443 1.648 1.852 2.034 2.158 2.145 1.818   0.302
[15,]   0.014 0.051 0.084 0.121 0.164 0.215 0.276 0.349 0.437 0.542 0.667 0.817 0.995 1.208 1.460 1.754 2.088 2.442 2.750 2.762   0.602
[16,]   0.012 0.033 0.050 0.070 0.093 0.121 0.156 0.198 0.251 0.317 0.400 0.507 0.644 0.822 1.059 1.377 1.815 2.425 3.283 4.375   1.817
[17,]   0.021 0.012 0.013 0.015 0.017 0.020 0.023 0.028 0.033 0.039 0.048 0.059 0.074 0.096 0.129 0.180 0.269 0.443 0.869 2.538 118.692
> ## tCopula -- df=4
> (theta <- iTau(tCopula(df=4), tau))
[1] 0.7071068
> dcop <- matrix(dCopula(umat, tCopula(param=theta, dim = 2, df=4)),
+                n1, n2)
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( tCopula(%.4g, df = 4) )", theta))
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( tCopula(%.4g, df = 4) )", theta))
> round(dcop, 3)
         [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]   [,21]
 [1,] 197.494 1.270 0.446 0.249 0.168 0.126 0.101 0.084 0.073 0.065 0.060 0.056 0.054 0.053 0.053 0.054 0.058 0.065 0.080 0.121   1.169
 [2,]   0.902 5.039 3.578 2.392 1.651 1.186 0.883 0.678 0.534 0.430 0.352 0.293 0.248 0.212 0.183 0.160 0.141 0.126 0.114 0.106   0.105
 [3,]   0.323 2.730 3.028 2.672 2.190 1.751 1.393 1.111 0.892 0.722 0.589 0.484 0.401 0.334 0.279 0.234 0.197 0.165 0.138 0.114   0.071
 [4,]   0.183 1.598 2.198 2.333 2.210 1.972 1.701 1.439 1.204 1.002 0.831 0.688 0.569 0.470 0.387 0.319 0.260 0.211 0.167 0.127   0.059
 [5,]   0.126 1.026 1.570 1.876 1.983 1.946 1.818 1.641 1.445 1.250 1.067 0.900 0.753 0.625 0.514 0.419 0.337 0.266 0.204 0.145   0.054
 [6,]   0.096 0.707 1.139 1.464 1.676 1.777 1.783 1.715 1.595 1.445 1.280 1.112 0.951 0.800 0.663 0.540 0.432 0.336 0.250 0.169   0.053
 [7,]   0.078 0.514 0.844 1.134 1.372 1.544 1.646 1.678 1.649 1.569 1.451 1.308 1.152 0.992 0.835 0.687 0.549 0.424 0.309 0.202   0.053
 [8,]   0.067 0.389 0.639 0.880 1.103 1.298 1.453 1.560 1.613 1.612 1.562 1.469 1.343 1.193 1.030 0.862 0.696 0.537 0.388 0.245   0.055
 [9,]   0.060 0.305 0.494 0.686 0.879 1.067 1.240 1.389 1.503 1.576 1.601 1.576 1.503 1.389 1.240 1.067 0.879 0.686 0.494 0.305   0.060
[10,]   0.055 0.245 0.388 0.537 0.696 0.862 1.030 1.193 1.343 1.469 1.562 1.612 1.613 1.560 1.453 1.298 1.103 0.880 0.639 0.389   0.067
[11,]   0.053 0.202 0.309 0.424 0.549 0.687 0.835 0.992 1.152 1.308 1.451 1.569 1.649 1.678 1.646 1.544 1.372 1.134 0.844 0.514   0.078
[12,]   0.053 0.169 0.250 0.336 0.432 0.540 0.663 0.800 0.951 1.112 1.280 1.445 1.595 1.715 1.783 1.777 1.676 1.464 1.139 0.707   0.096
[13,]   0.054 0.145 0.204 0.266 0.337 0.419 0.514 0.625 0.753 0.900 1.067 1.250 1.445 1.641 1.818 1.946 1.983 1.876 1.570 1.026   0.126
[14,]   0.059 0.127 0.167 0.211 0.260 0.319 0.387 0.470 0.569 0.688 0.831 1.002 1.204 1.439 1.701 1.972 2.210 2.333 2.198 1.598   0.183
[15,]   0.071 0.114 0.138 0.165 0.197 0.234 0.279 0.334 0.401 0.484 0.589 0.722 0.892 1.111 1.393 1.751 2.190 2.672 3.028 2.730   0.323
[16,]   0.105 0.106 0.114 0.126 0.141 0.160 0.183 0.212 0.248 0.293 0.352 0.430 0.534 0.678 0.883 1.186 1.651 2.392 3.578 5.039   0.902
[17,]   1.169 0.121 0.080 0.065 0.058 0.054 0.053 0.053 0.054 0.056 0.060 0.065 0.073 0.084 0.101 0.126 0.168 0.249 0.446 1.270 197.494
> 
> ## galambosCopula
> #
> (theta <- iTau(galambosCopula(), tau))
[1] 1.284823
> stopifnot(all.equal(tau, tau(galambosCopula(theta)), tolerance = 1e-5))
> dcop <- matrix(dCopula(umat, galambosCopula(param=theta)),
+                n1, n2)
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( galambosCopula(%.4g) )", theta))
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( galambosCopula(%.4g) )", theta))
> round(dcop, 3)
        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]   [,21]
 [1,] 30.876 4.566 2.428 1.545 1.066 0.770 0.572 0.433 0.331 0.255 0.197 0.151 0.115 0.087 0.064 0.046 0.031 0.020 0.011 0.004   0.000
 [2,]  3.787 3.300 2.698 2.232 1.858 1.550 1.293 1.075 0.891 0.734 0.599 0.484 0.385 0.301 0.230 0.169 0.118 0.076 0.042 0.016   0.000
 [3,]  1.910 2.479 2.326 2.120 1.904 1.692 1.489 1.298 1.121 0.957 0.808 0.671 0.549 0.439 0.341 0.256 0.182 0.118 0.066 0.026   0.000
 [4,]  1.164 1.914 1.972 1.927 1.834 1.713 1.576 1.430 1.280 1.129 0.982 0.839 0.703 0.575 0.457 0.348 0.251 0.166 0.094 0.037   0.000
 [5,]  0.770 1.495 1.655 1.710 1.708 1.665 1.594 1.500 1.389 1.266 1.133 0.996 0.856 0.717 0.582 0.453 0.333 0.224 0.129 0.051   0.000
 [6,]  0.533 1.172 1.375 1.489 1.551 1.572 1.561 1.520 1.455 1.368 1.263 1.142 1.009 0.867 0.721 0.574 0.430 0.295 0.172 0.069   0.000
 [7,]  0.378 0.917 1.129 1.274 1.377 1.446 1.485 1.496 1.479 1.435 1.366 1.273 1.158 1.024 0.875 0.715 0.549 0.384 0.229 0.093   0.001
 [8,]  0.272 0.713 0.915 1.070 1.196 1.297 1.376 1.431 1.460 1.463 1.438 1.384 1.300 1.186 1.045 0.880 0.695 0.500 0.305 0.126   0.001
 [9,]  0.197 0.549 0.730 0.880 1.014 1.133 1.239 1.328 1.399 1.448 1.472 1.465 1.424 1.346 1.228 1.071 0.875 0.650 0.409 0.174   0.001
[10,]  0.141 0.417 0.570 0.707 0.836 0.961 1.080 1.193 1.297 1.387 1.458 1.505 1.519 1.493 1.418 1.288 1.098 0.849 0.555 0.244   0.002
[11,]  0.100 0.309 0.433 0.550 0.667 0.786 0.908 1.032 1.155 1.277 1.390 1.490 1.566 1.608 1.600 1.526 1.369 1.115 0.766 0.353   0.002
[12,]  0.069 0.221 0.318 0.412 0.510 0.615 0.728 0.850 0.980 1.118 1.262 1.406 1.544 1.664 1.746 1.766 1.687 1.471 1.083 0.533   0.004
[13,]  0.046 0.151 0.221 0.292 0.368 0.453 0.548 0.656 0.778 0.916 1.071 1.243 1.429 1.623 1.810 1.963 2.031 1.936 1.571 0.855   0.007
[14,]  0.028 0.095 0.141 0.189 0.243 0.304 0.376 0.460 0.560 0.679 0.822 0.993 1.198 1.441 1.721 2.026 2.316 2.491 2.333 1.489   0.013
[15,]  0.015 0.052 0.078 0.105 0.137 0.174 0.219 0.273 0.340 0.424 0.530 0.667 0.844 1.076 1.385 1.794 2.323 2.946 3.436 2.914   0.035
[16,]  0.006 0.019 0.029 0.040 0.053 0.068 0.087 0.110 0.140 0.179 0.230 0.299 0.396 0.535 0.742 1.067 1.600 2.525 4.166 6.313   0.174
[17,]  0.000 0.000 0.000 0.000 0.000 0.000 0.000 0.001 0.001 0.001 0.001 0.002 0.002 0.003 0.004 0.007 0.012 0.023 0.059 0.290 333.569
> 
> ## plackettCopula
> #
> (theta <- iTau(plackettCopula(), tau))
[1] 11.39548
> stopifnot(all.equal(tau, tau(plackettCopula(theta)), tolerance = 1e-5))
> dcop <- matrix(dCopula(umat, plackettCopula(param=theta)),
+                n1, n2)
> filled.contour(u1, u2, log(dcop), color.palette = fCols,
+                main=sprintf("log Density( plackettCopula(%.4g) )", theta))
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( plackettCopula(%.4g) )", theta))
> round(dcop, 3)
        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20]  [,21]
 [1,] 10.943 4.938 2.755 1.752 1.210 0.886 0.676 0.533 0.430 0.355 0.298 0.254 0.218 0.190 0.167 0.148 0.132 0.118 0.106 0.097  0.088
 [2,]  4.199 4.139 3.305 2.424 1.753 1.287 0.968 0.746 0.589 0.474 0.389 0.324 0.273 0.234 0.202 0.176 0.155 0.137 0.122 0.109  0.099
 [3,]  2.169 2.752 2.937 2.686 2.218 1.739 1.340 1.035 0.808 0.640 0.515 0.421 0.349 0.292 0.248 0.213 0.184 0.161 0.142 0.125  0.112
 [4,]  1.319 1.784 2.203 2.415 2.351 2.078 1.723 1.380 1.091 0.862 0.687 0.553 0.450 0.371 0.310 0.261 0.223 0.191 0.166 0.145  0.128
 [5,]  0.886 1.199 1.561 1.899 2.115 2.139 1.980 1.712 1.415 1.142 0.913 0.731 0.588 0.478 0.392 0.325 0.273 0.231 0.197 0.170  0.148
 [6,]  0.635 0.843 1.108 1.409 1.702 1.920 1.996 1.913 1.710 1.451 1.192 0.963 0.774 0.623 0.504 0.412 0.339 0.282 0.237 0.201  0.172
 [7,]  0.477 0.619 0.804 1.032 1.295 1.563 1.782 1.895 1.869 1.718 1.492 1.246 1.015 0.818 0.657 0.530 0.430 0.351 0.290 0.241  0.203
 [8,]  0.372 0.470 0.600 0.765 0.969 1.206 1.457 1.681 1.823 1.844 1.738 1.542 1.305 1.071 0.864 0.692 0.554 0.446 0.361 0.295  0.244
 [9,]  0.298 0.368 0.460 0.578 0.728 0.913 1.133 1.374 1.603 1.773 1.836 1.773 1.603 1.374 1.133 0.913 0.728 0.578 0.460 0.368  0.298
[10,]  0.244 0.295 0.361 0.446 0.554 0.692 0.864 1.071 1.305 1.542 1.738 1.844 1.823 1.681 1.457 1.206 0.969 0.765 0.600 0.470  0.372
[11,]  0.203 0.241 0.290 0.351 0.430 0.530 0.657 0.818 1.015 1.246 1.492 1.718 1.869 1.895 1.782 1.563 1.295 1.032 0.804 0.619  0.477
[12,]  0.172 0.201 0.237 0.282 0.339 0.412 0.504 0.623 0.774 0.963 1.192 1.451 1.710 1.913 1.996 1.920 1.702 1.409 1.108 0.843  0.635
[13,]  0.148 0.170 0.197 0.231 0.273 0.325 0.392 0.478 0.588 0.731 0.913 1.142 1.415 1.712 1.980 2.139 2.115 1.899 1.561 1.199  0.886
[14,]  0.128 0.145 0.166 0.191 0.223 0.261 0.310 0.371 0.450 0.553 0.687 0.862 1.091 1.380 1.723 2.078 2.351 2.415 2.203 1.784  1.319
[15,]  0.112 0.125 0.142 0.161 0.184 0.213 0.248 0.292 0.349 0.421 0.515 0.640 0.808 1.035 1.340 1.739 2.218 2.686 2.937 2.752  2.169
[16,]  0.099 0.109 0.122 0.137 0.155 0.176 0.202 0.234 0.273 0.324 0.389 0.474 0.589 0.746 0.968 1.287 1.753 2.424 3.305 4.139  4.199
[17,]  0.088 0.097 0.106 0.118 0.132 0.148 0.167 0.190 0.218 0.254 0.298 0.355 0.430 0.533 0.676 0.886 1.210 1.752 2.755 4.938 10.943
> 
> 
> ## amhCopula
> tau <- 0.3 ## -- to be in range for AMH
> (theta <- iTau(amhCopula(), tau))
[1] 0.9429734
> stopifnot(all.equal(tau, tau(amhCopula(theta)), tolerance = 1e-5))
> dcop <- matrix(dCopula(umat, amhCopula(param=theta)),
+                n1, n2)
> filled.contour(u1, u2, dcop, color.palette = fCols,
+                main=sprintf("Density( amhCopula(%.4g) )", theta))
> round(dcop, 3)
        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 16.441 5.247 2.516 1.472 0.965 0.681 0.506 0.391 0.311 0.253 0.210 0.177 0.152 0.131 0.114 0.101 0.089 0.080 0.072 0.065 0.059
 [2,]  4.251 3.657 2.746 2.074 1.605 1.272 1.031 0.851 0.714 0.607 0.522 0.454 0.398 0.352 0.314 0.281 0.253 0.229 0.209 0.191 0.175
 [3,]  1.890 2.388 2.243 1.967 1.694 1.457 1.258 1.093 0.956 0.842 0.747 0.666 0.597 0.539 0.488 0.444 0.406 0.372 0.342 0.316 0.293
 [4,]  1.064 1.647 1.774 1.719 1.600 1.464 1.329 1.205 1.092 0.992 0.903 0.824 0.755 0.693 0.638 0.589 0.545 0.506 0.471 0.439 0.411
 [5,]  0.681 1.197 1.413 1.473 1.454 1.397 1.322 1.242 1.161 1.083 1.009 0.941 0.878 0.820 0.767 0.718 0.673 0.632 0.595 0.560 0.529
 [6,]  0.473 0.906 1.144 1.260 1.303 1.302 1.276 1.235 1.186 1.133 1.079 1.026 0.974 0.924 0.877 0.833 0.790 0.751 0.714 0.679 0.647
 [7,]  0.347 0.709 0.942 1.083 1.162 1.201 1.212 1.204 1.184 1.156 1.122 1.086 1.048 1.010 0.972 0.934 0.898 0.862 0.828 0.795 0.765
 [8,]  0.266 0.570 0.787 0.936 1.037 1.102 1.140 1.159 1.164 1.159 1.146 1.127 1.104 1.079 1.052 1.024 0.995 0.967 0.938 0.910 0.883
 [9,]  0.210 0.467 0.666 0.816 0.928 1.009 1.068 1.108 1.134 1.148 1.154 1.153 1.146 1.135 1.120 1.104 1.085 1.065 1.044 1.022 1.000
[10,]  0.170 0.390 0.571 0.716 0.833 0.925 0.997 1.054 1.097 1.129 1.151 1.166 1.175 1.179 1.178 1.174 1.167 1.157 1.145 1.132 1.118
[11,]  0.141 0.331 0.495 0.633 0.750 0.848 0.931 0.999 1.056 1.103 1.141 1.171 1.195 1.213 1.226 1.236 1.241 1.243 1.243 1.240 1.236
[12,]  0.118 0.284 0.432 0.563 0.679 0.780 0.868 0.946 1.014 1.073 1.124 1.169 1.207 1.239 1.267 1.290 1.309 1.325 1.337 1.347 1.353
[13,]  0.101 0.246 0.381 0.504 0.616 0.718 0.811 0.895 0.971 1.041 1.104 1.160 1.212 1.258 1.300 1.337 1.371 1.401 1.427 1.451 1.471
[14,]  0.087 0.215 0.338 0.454 0.562 0.663 0.758 0.846 0.929 1.007 1.080 1.148 1.212 1.271 1.327 1.379 1.427 1.472 1.514 1.553 1.589
[15,]  0.076 0.190 0.302 0.410 0.514 0.613 0.709 0.800 0.889 0.973 1.054 1.132 1.207 1.279 1.348 1.414 1.478 1.539 1.597 1.653 1.706
[16,]  0.067 0.169 0.272 0.373 0.472 0.569 0.664 0.757 0.849 0.939 1.028 1.114 1.199 1.283 1.365 1.445 1.524 1.601 1.677 1.752 1.824
[17,]  0.059 0.152 0.246 0.340 0.435 0.529 0.623 0.718 0.812 0.906 1.000 1.095 1.189 1.283 1.377 1.471 1.565 1.659 1.753 1.847 1.939
> 
> showProc.time()
Time (user system elapsed): 0.14 0.02 0.15 
> if(!dev.interactive()) dev.off()
null device 
          1 
> 
> 
> ### d > 2 ######################################################################
> 
> ### for package fCopulae
> ## dcop.f <- darchmCopula(umat, alpha=2, type=1, alternative=T)
> ## round(matrix(dcop.f, n1, n2), 3)
> ## dcop.fCopulae <- darchmCopula(umat, alpha=theta.fr, type = 5, alternative = TRUE)
> ## round (matrix(dcop.fCopulae, n1, n2), 3)
> 
> ## dim = 3, frankCopula
> um3 <- Exp.grid(u1=u1, u2=u2, u3=u1)
> dcE <- tryCatch(dcopula(um3), error=identity)
> dcop <- dCopula(um3, frankCopula(param=theta.fr, dim = 3))
> stopifnot(dcop >= 0, # no NA here
+           inherits(dcE, "error"),
+           grepl("defunct", conditionMessage(dcE), ignore.case=TRUE))
> round(array(dcop, c(n1,n2,n1)), 3)
, , 1

        [,1]   [,2]   [,3]   [,4]   [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 32.553 24.578 18.450 13.849 10.396 7.804 5.858 4.397 3.301 2.478 1.860 1.396 1.048 0.787 0.591 0.443 0.333 0.250 0.188 0.141 0.106
 [2,] 22.877 17.301 13.004  9.771  7.340 5.513 4.140 3.109 2.334 1.752 1.316 0.988 0.741 0.557 0.418 0.314 0.235 0.177 0.133 0.100 0.075
 [3,] 15.985 12.103  9.105  6.846  5.145 3.866 2.904 2.181 1.638 1.230 0.923 0.693 0.520 0.391 0.293 0.220 0.165 0.124 0.093 0.070 0.053
 [4,] 11.169  8.464  6.371  4.793  3.603 2.708 2.035 1.528 1.148 0.862 0.647 0.486 0.365 0.274 0.206 0.154 0.116 0.087 0.065 0.049 0.037
 [5,]  7.804  5.917  4.456  3.353  2.522 1.896 1.424 1.070 0.804 0.604 0.453 0.340 0.256 0.192 0.144 0.108 0.081 0.061 0.046 0.034 0.026
 [6,]  5.453  4.136  3.116  2.345  1.764 1.326 0.997 0.749 0.562 0.422 0.317 0.238 0.179 0.134 0.101 0.076 0.057 0.043 0.032 0.024 0.018
 [7,]  3.810  2.891  2.178  1.640  1.234 0.927 0.697 0.524 0.393 0.295 0.222 0.167 0.125 0.094 0.070 0.053 0.040 0.030 0.022 0.017 0.013
 [8,]  2.662  2.020  1.522  1.146  0.862 0.648 0.487 0.366 0.275 0.207 0.155 0.116 0.087 0.066 0.049 0.037 0.028 0.021 0.016 0.012 0.009
 [9,]  1.860  1.412  1.064  0.801  0.603 0.453 0.341 0.256 0.192 0.144 0.108 0.081 0.061 0.046 0.034 0.026 0.019 0.015 0.011 0.008 0.006
[10,]  1.300  0.986  0.744  0.560  0.421 0.317 0.238 0.179 0.134 0.101 0.076 0.057 0.043 0.032 0.024 0.018 0.014 0.010 0.008 0.006 0.004
[11,]  0.908  0.689  0.520  0.391  0.294 0.221 0.166 0.125 0.094 0.071 0.053 0.040 0.030 0.022 0.017 0.013 0.009 0.007 0.005 0.004 0.003
[12,]  0.634  0.482  0.363  0.273  0.206 0.155 0.116 0.087 0.066 0.049 0.037 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004 0.003 0.002
[13,]  0.443  0.337  0.254  0.191  0.144 0.108 0.081 0.061 0.046 0.034 0.026 0.019 0.015 0.011 0.008 0.006 0.005 0.003 0.003 0.002 0.001
[14,]  0.310  0.235  0.177  0.134  0.100 0.076 0.057 0.043 0.032 0.024 0.018 0.014 0.010 0.008 0.006 0.004 0.003 0.002 0.002 0.001 0.001
[15,]  0.216  0.164  0.124  0.093  0.070 0.053 0.040 0.030 0.022 0.017 0.013 0.009 0.007 0.005 0.004 0.003 0.002 0.002 0.001 0.001 0.001
[16,]  0.151  0.115  0.087  0.065  0.049 0.037 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004 0.003 0.002 0.002 0.001 0.001 0.001 0.001
[17,]  0.106  0.081  0.061  0.046  0.034 0.026 0.019 0.015 0.011 0.008 0.006 0.005 0.003 0.003 0.002 0.001 0.001 0.001 0.001 0.000 0.000

, , 2

        [,1]   [,2]   [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 22.877 17.301 13.004 9.771 7.340 5.513 4.140 3.109 2.334 1.752 1.316 0.988 0.741 0.557 0.418 0.314 0.235 0.177 0.133 0.100 0.075
 [2,] 16.109 13.302 10.704 8.469 6.616 5.117 3.929 3.000 2.281 1.729 1.307 0.987 0.744 0.560 0.421 0.317 0.238 0.179 0.134 0.101 0.076
 [3,] 11.272  9.915  8.389 6.899 5.550 4.391 3.430 2.653 2.037 1.556 1.183 0.896 0.678 0.512 0.386 0.290 0.218 0.164 0.123 0.093 0.070
 [4,]  7.884  7.251  6.358 5.377 4.421 3.557 2.814 2.198 1.700 1.306 0.997 0.758 0.575 0.434 0.328 0.247 0.186 0.140 0.105 0.079 0.060
 [5,]  5.513  5.231  4.706 4.060 3.391 2.762 2.206 1.736 1.350 1.041 0.797 0.608 0.461 0.349 0.264 0.199 0.150 0.113 0.085 0.064 0.048
 [6,]  3.854  3.738  3.424 2.997 2.532 2.080 1.673 1.323 1.033 0.799 0.614 0.469 0.356 0.270 0.204 0.154 0.116 0.087 0.066 0.049 0.037
 [7,]  2.694  2.653  2.461 2.176 1.853 1.533 1.238 0.983 0.770 0.597 0.459 0.351 0.267 0.203 0.153 0.116 0.087 0.066 0.049 0.037 0.028
 [8,]  1.882  1.874  1.754 1.563 1.338 1.112 0.901 0.718 0.563 0.437 0.337 0.258 0.196 0.149 0.113 0.085 0.064 0.048 0.036 0.027 0.021
 [9,]  1.316  1.320  1.243 1.113 0.957 0.797 0.648 0.517 0.406 0.316 0.243 0.186 0.142 0.108 0.082 0.062 0.047 0.035 0.026 0.020 0.015
[10,]  0.919  0.927  0.877 0.788 0.680 0.568 0.462 0.369 0.290 0.226 0.174 0.133 0.102 0.077 0.059 0.044 0.033 0.025 0.019 0.014 0.011
[11,]  0.642  0.650  0.617 0.556 0.480 0.402 0.328 0.262 0.206 0.161 0.124 0.095 0.072 0.055 0.042 0.031 0.024 0.018 0.013 0.010 0.008
[12,]  0.449  0.456  0.433 0.391 0.338 0.283 0.231 0.185 0.146 0.113 0.088 0.067 0.051 0.039 0.029 0.022 0.017 0.013 0.010 0.007 0.005
[13,]  0.314  0.319  0.304 0.274 0.238 0.199 0.163 0.130 0.103 0.080 0.062 0.047 0.036 0.027 0.021 0.016 0.012 0.009 0.007 0.005 0.004
[14,]  0.219  0.223  0.213 0.192 0.167 0.140 0.114 0.091 0.072 0.056 0.043 0.033 0.025 0.019 0.015 0.011 0.008 0.006 0.005 0.004 0.003
[15,]  0.153  0.156  0.149 0.135 0.117 0.098 0.080 0.064 0.051 0.039 0.030 0.023 0.018 0.014 0.010 0.008 0.006 0.004 0.003 0.002 0.002
[16,]  0.107  0.109  0.104 0.094 0.082 0.069 0.056 0.045 0.035 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004 0.003 0.002 0.002 0.001
[17,]  0.075  0.077  0.073 0.066 0.058 0.048 0.039 0.032 0.025 0.019 0.015 0.012 0.009 0.007 0.005 0.004 0.003 0.002 0.002 0.001 0.001

, , 3

        [,1]   [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 15.985 12.103 9.105 6.846 5.145 3.866 2.904 2.181 1.638 1.230 0.923 0.693 0.520 0.391 0.293 0.220 0.165 0.124 0.093 0.070 0.053
 [2,] 11.272  9.915 8.389 6.899 5.550 4.391 3.430 2.653 2.037 1.556 1.183 0.896 0.678 0.512 0.386 0.290 0.218 0.164 0.123 0.093 0.070
 [3,]  7.896  7.744 7.157 6.311 5.362 4.426 3.571 2.832 2.217 1.717 1.320 1.009 0.768 0.582 0.440 0.332 0.251 0.189 0.142 0.107 0.081
 [4,]  5.526  5.854 5.772 5.365 4.755 4.058 3.362 2.721 2.163 1.696 1.316 1.013 0.775 0.590 0.448 0.339 0.256 0.193 0.145 0.109 0.083
 [5,]  3.866  4.325 4.468 4.318 3.951 3.461 2.927 2.408 1.939 1.535 1.200 0.929 0.713 0.545 0.414 0.314 0.238 0.179 0.135 0.102 0.077
 [6,]  2.704  3.143 3.358 3.338 3.128 2.794 2.400 2.000 1.626 1.297 1.020 0.793 0.611 0.468 0.357 0.271 0.205 0.155 0.117 0.088 0.067
 [7,]  1.890  2.257 2.470 2.506 2.390 2.166 1.883 1.583 1.297 1.041 0.822 0.641 0.496 0.380 0.290 0.221 0.167 0.126 0.095 0.072 0.054
 [8,]  1.321  1.608 1.789 1.843 1.779 1.629 1.429 1.210 0.996 0.803 0.636 0.498 0.385 0.296 0.226 0.172 0.131 0.099 0.075 0.056 0.043
 [9,]  0.923  1.139 1.282 1.335 1.300 1.200 1.059 0.901 0.745 0.602 0.479 0.375 0.291 0.224 0.171 0.130 0.099 0.075 0.056 0.043 0.032
[10,]  0.645  0.803 0.912 0.956 0.938 0.870 0.771 0.659 0.546 0.443 0.352 0.276 0.215 0.165 0.126 0.096 0.073 0.055 0.042 0.032 0.024
[11,]  0.451  0.565 0.645 0.680 0.670 0.624 0.555 0.475 0.395 0.320 0.255 0.201 0.156 0.120 0.092 0.070 0.053 0.040 0.030 0.023 0.017
[12,]  0.315  0.397 0.455 0.481 0.475 0.444 0.396 0.339 0.282 0.230 0.183 0.144 0.112 0.086 0.066 0.050 0.038 0.029 0.022 0.016 0.012
[13,]  0.220  0.278 0.320 0.339 0.336 0.314 0.280 0.241 0.201 0.163 0.130 0.102 0.080 0.061 0.047 0.036 0.027 0.021 0.016 0.012 0.009
[14,]  0.154  0.195 0.224 0.238 0.236 0.222 0.198 0.170 0.142 0.115 0.092 0.073 0.056 0.044 0.033 0.025 0.019 0.015 0.011 0.008 0.006
[15,]  0.108  0.136 0.157 0.167 0.166 0.156 0.139 0.120 0.100 0.081 0.065 0.051 0.040 0.031 0.024 0.018 0.014 0.010 0.008 0.006 0.004
[16,]  0.075  0.095 0.110 0.117 0.117 0.109 0.098 0.084 0.070 0.057 0.046 0.036 0.028 0.022 0.017 0.013 0.010 0.007 0.005 0.004 0.003
[17,]  0.053  0.067 0.077 0.083 0.082 0.077 0.069 0.059 0.050 0.040 0.032 0.025 0.020 0.015 0.012 0.009 0.007 0.005 0.004 0.003 0.002

, , 4

        [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 11.169 8.464 6.371 4.793 3.603 2.708 2.035 1.528 1.148 0.862 0.647 0.486 0.365 0.274 0.206 0.154 0.116 0.087 0.065 0.049 0.037
 [2,]  7.884 7.251 6.358 5.377 4.421 3.557 2.814 2.198 1.700 1.306 0.997 0.758 0.575 0.434 0.328 0.247 0.186 0.140 0.105 0.079 0.060
 [3,]  5.526 5.854 5.772 5.365 4.755 4.058 3.362 2.721 2.163 1.696 1.316 1.013 0.775 0.590 0.448 0.339 0.256 0.193 0.145 0.109 0.083
 [4,]  3.870 4.534 4.878 4.879 4.596 4.122 3.553 2.968 2.418 1.933 1.522 1.184 0.913 0.700 0.534 0.405 0.307 0.232 0.175 0.132 0.100
 [5,]  2.708 3.408 3.909 4.135 4.086 3.815 3.399 2.917 2.427 1.973 1.573 1.237 0.961 0.741 0.567 0.432 0.328 0.249 0.188 0.142 0.107
 [6,]  1.894 2.507 3.013 3.322 3.406 3.282 3.004 2.634 2.231 1.838 1.482 1.175 0.920 0.712 0.548 0.418 0.318 0.241 0.183 0.138 0.105
 [7,]  1.325 1.816 2.257 2.566 2.704 2.670 2.494 2.225 1.911 1.593 1.296 1.034 0.814 0.633 0.488 0.374 0.285 0.216 0.164 0.124 0.094
 [8,]  0.926 1.302 1.656 1.926 2.071 2.082 1.975 1.785 1.550 1.303 1.067 0.857 0.677 0.528 0.408 0.313 0.239 0.182 0.138 0.104 0.079
 [9,]  0.647 0.926 1.198 1.415 1.544 1.573 1.510 1.378 1.207 1.021 0.841 0.678 0.537 0.420 0.325 0.250 0.191 0.145 0.110 0.083 0.063
[10,]  0.452 0.655 0.858 1.025 1.130 1.162 1.125 1.034 0.911 0.774 0.640 0.518 0.411 0.322 0.250 0.192 0.147 0.112 0.085 0.064 0.049
[11,]  0.316 0.462 0.610 0.734 0.816 0.845 0.823 0.760 0.672 0.574 0.476 0.385 0.307 0.241 0.187 0.144 0.110 0.084 0.064 0.048 0.037
[12,]  0.221 0.325 0.431 0.522 0.583 0.607 0.594 0.550 0.488 0.418 0.347 0.282 0.224 0.176 0.137 0.106 0.081 0.062 0.047 0.035 0.027
[13,]  0.154 0.228 0.304 0.369 0.414 0.432 0.424 0.394 0.351 0.301 0.250 0.203 0.162 0.127 0.099 0.076 0.058 0.045 0.034 0.026 0.019
[14,]  0.108 0.160 0.214 0.260 0.293 0.306 0.301 0.281 0.250 0.214 0.179 0.145 0.116 0.091 0.071 0.055 0.042 0.032 0.024 0.018 0.014
[15,]  0.075 0.112 0.150 0.183 0.206 0.216 0.213 0.198 0.177 0.152 0.127 0.103 0.082 0.065 0.050 0.039 0.030 0.023 0.017 0.013 0.010
[16,]  0.053 0.078 0.105 0.128 0.145 0.152 0.150 0.140 0.125 0.107 0.089 0.073 0.058 0.046 0.036 0.027 0.021 0.016 0.012 0.009 0.007
[17,]  0.037 0.055 0.074 0.091 0.102 0.107 0.106 0.099 0.088 0.076 0.063 0.052 0.041 0.032 0.025 0.019 0.015 0.011 0.009 0.007 0.005

, , 5

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 7.804 5.917 4.456 3.353 2.522 1.896 1.424 1.070 0.804 0.604 0.453 0.340 0.256 0.192 0.144 0.108 0.081 0.061 0.046 0.034 0.026
 [2,] 5.513 5.231 4.706 4.060 3.391 2.762 2.206 1.736 1.350 1.041 0.797 0.608 0.461 0.349 0.264 0.199 0.150 0.113 0.085 0.064 0.048
 [3,] 3.866 4.325 4.468 4.318 3.951 3.461 2.927 2.408 1.939 1.535 1.200 0.929 0.713 0.545 0.414 0.314 0.238 0.179 0.135 0.102 0.077
 [4,] 2.708 3.408 3.909 4.135 4.086 3.815 3.399 2.917 2.427 1.973 1.573 1.237 0.961 0.741 0.567 0.432 0.328 0.249 0.188 0.142 0.107
 [5,] 1.896 2.594 3.216 3.649 3.839 3.784 3.530 3.145 2.699 2.248 1.827 1.458 1.147 0.892 0.687 0.527 0.401 0.305 0.231 0.174 0.132
 [6,] 1.326 1.925 2.526 3.024 3.342 3.443 3.338 3.072 2.709 2.306 1.907 1.543 1.226 0.961 0.746 0.574 0.439 0.334 0.254 0.192 0.146
 [7,] 0.927 1.404 1.919 2.391 2.743 2.926 2.926 2.766 2.495 2.163 1.816 1.487 1.193 0.942 0.735 0.568 0.436 0.333 0.253 0.192 0.146
 [8,] 0.648 1.011 1.423 1.825 2.154 2.359 2.416 2.333 2.142 1.886 1.603 1.325 1.072 0.851 0.667 0.518 0.398 0.305 0.232 0.176 0.134
 [9,] 0.453 0.721 1.037 1.358 1.636 1.827 1.906 1.871 1.742 1.552 1.332 1.110 0.903 0.721 0.567 0.441 0.340 0.261 0.199 0.151 0.115
[10,] 0.317 0.511 0.746 0.992 1.213 1.375 1.453 1.444 1.359 1.221 1.056 0.885 0.724 0.580 0.458 0.357 0.276 0.211 0.161 0.123 0.093
[11,] 0.221 0.361 0.532 0.715 0.884 1.012 1.081 1.083 1.027 0.930 0.808 0.681 0.558 0.449 0.355 0.277 0.214 0.165 0.126 0.096 0.073
[12,] 0.155 0.254 0.377 0.511 0.636 0.734 0.789 0.796 0.759 0.690 0.603 0.509 0.419 0.337 0.267 0.209 0.162 0.124 0.095 0.072 0.055
[13,] 0.108 0.178 0.266 0.363 0.454 0.527 0.569 0.577 0.552 0.504 0.441 0.374 0.308 0.248 0.197 0.154 0.119 0.092 0.070 0.053 0.041
[14,] 0.076 0.125 0.187 0.256 0.322 0.375 0.406 0.413 0.397 0.363 0.319 0.270 0.223 0.180 0.143 0.112 0.087 0.067 0.051 0.039 0.030
[15,] 0.053 0.088 0.132 0.180 0.227 0.265 0.288 0.294 0.283 0.259 0.228 0.193 0.160 0.129 0.103 0.080 0.062 0.048 0.037 0.028 0.021
[16,] 0.037 0.061 0.092 0.127 0.160 0.187 0.204 0.208 0.200 0.184 0.162 0.137 0.114 0.092 0.073 0.057 0.044 0.034 0.026 0.020 0.015
[17,] 0.026 0.043 0.065 0.089 0.113 0.132 0.144 0.147 0.142 0.131 0.115 0.098 0.081 0.065 0.052 0.041 0.032 0.024 0.019 0.014 0.011

, , 6

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 5.453 4.136 3.116 2.345 1.764 1.326 0.997 0.749 0.562 0.422 0.317 0.238 0.179 0.134 0.101 0.076 0.057 0.043 0.032 0.024 0.018
 [2,] 3.854 3.738 3.424 2.997 2.532 2.080 1.673 1.323 1.033 0.799 0.614 0.469 0.356 0.270 0.204 0.154 0.116 0.087 0.066 0.049 0.037
 [3,] 2.704 3.143 3.358 3.338 3.128 2.794 2.400 2.000 1.626 1.297 1.020 0.793 0.611 0.468 0.357 0.271 0.205 0.155 0.117 0.088 0.067
 [4,] 1.894 2.507 3.013 3.322 3.406 3.282 3.004 2.634 2.231 1.838 1.482 1.175 0.920 0.712 0.548 0.418 0.318 0.241 0.183 0.138 0.105
 [5,] 1.326 1.925 2.526 3.024 3.342 3.443 3.338 3.072 2.709 2.306 1.907 1.543 1.226 0.961 0.746 0.574 0.439 0.334 0.254 0.192 0.146
 [6,] 0.928 1.438 2.014 2.568 3.014 3.284 3.348 3.220 2.946 2.585 2.192 1.809 1.460 1.159 0.907 0.703 0.541 0.414 0.315 0.239 0.182
 [7,] 0.649 1.054 1.546 2.068 2.544 2.899 3.081 3.076 2.908 2.623 2.275 1.911 1.565 1.256 0.993 0.775 0.599 0.460 0.351 0.267 0.203
 [8,] 0.454 0.761 1.156 1.600 2.040 2.407 2.645 2.723 2.644 2.441 2.159 1.843 1.529 1.240 0.987 0.775 0.602 0.464 0.355 0.271 0.207
 [9,] 0.317 0.544 0.847 1.203 1.574 1.907 2.150 2.267 2.250 2.117 1.903 1.646 1.380 1.129 0.905 0.714 0.557 0.430 0.330 0.252 0.193
[10,] 0.222 0.387 0.612 0.885 1.180 1.459 1.677 1.801 1.817 1.735 1.579 1.381 1.168 0.962 0.775 0.614 0.481 0.372 0.286 0.219 0.167
[11,] 0.155 0.273 0.438 0.641 0.867 1.087 1.268 1.380 1.410 1.362 1.252 1.103 0.939 0.778 0.629 0.501 0.393 0.305 0.235 0.180 0.138
[12,] 0.108 0.192 0.311 0.460 0.628 0.795 0.937 1.030 1.063 1.035 0.958 0.849 0.727 0.604 0.491 0.391 0.307 0.239 0.184 0.141 0.108
[13,] 0.076 0.135 0.220 0.327 0.450 0.574 0.681 0.754 0.783 0.767 0.714 0.636 0.546 0.456 0.371 0.296 0.233 0.182 0.140 0.107 0.082
[14,] 0.053 0.095 0.155 0.231 0.320 0.410 0.490 0.545 0.569 0.559 0.523 0.467 0.402 0.336 0.274 0.219 0.173 0.135 0.104 0.080 0.061
[15,] 0.037 0.066 0.109 0.163 0.226 0.291 0.349 0.390 0.408 0.403 0.377 0.338 0.292 0.244 0.199 0.160 0.126 0.098 0.076 0.058 0.045
[16,] 0.026 0.047 0.076 0.115 0.160 0.206 0.247 0.277 0.291 0.287 0.270 0.242 0.209 0.175 0.143 0.115 0.091 0.071 0.055 0.042 0.032
[17,] 0.018 0.033 0.054 0.081 0.113 0.146 0.175 0.197 0.207 0.205 0.193 0.173 0.150 0.125 0.103 0.082 0.065 0.051 0.039 0.030 0.023

, , 7

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 3.810 2.891 2.178 1.640 1.234 0.927 0.697 0.524 0.393 0.295 0.222 0.167 0.125 0.094 0.070 0.053 0.040 0.030 0.022 0.017 0.013
 [2,] 2.694 2.653 2.461 2.176 1.853 1.533 1.238 0.983 0.770 0.597 0.459 0.351 0.267 0.203 0.153 0.116 0.087 0.066 0.049 0.037 0.028
 [3,] 1.890 2.257 2.470 2.506 2.390 2.166 1.883 1.583 1.297 1.041 0.822 0.641 0.496 0.380 0.290 0.221 0.167 0.126 0.095 0.072 0.054
 [4,] 1.325 1.816 2.257 2.566 2.704 2.670 2.494 2.225 1.911 1.593 1.296 1.034 0.814 0.633 0.488 0.374 0.285 0.216 0.164 0.124 0.094
 [5,] 0.927 1.404 1.919 2.391 2.743 2.926 2.926 2.766 2.495 2.163 1.816 1.487 1.193 0.942 0.735 0.568 0.436 0.333 0.253 0.192 0.146
 [6,] 0.649 1.054 1.546 2.068 2.544 2.899 3.081 3.076 2.908 2.623 2.275 1.911 1.565 1.256 0.993 0.775 0.599 0.460 0.351 0.267 0.203
 [7,] 0.454 0.774 1.197 1.689 2.195 2.641 2.957 3.098 3.056 2.861 2.560 2.207 1.845 1.505 1.204 0.949 0.740 0.571 0.438 0.334 0.255
 [8,] 0.317 0.561 0.899 1.321 1.790 2.249 2.627 2.866 2.934 2.838 2.612 2.305 1.965 1.628 1.319 1.049 0.823 0.639 0.492 0.377 0.289
 [9,] 0.222 0.402 0.662 1.000 1.399 1.816 2.194 2.472 2.610 2.595 2.447 2.204 1.911 1.605 1.314 1.055 0.834 0.651 0.503 0.386 0.296
[10,] 0.155 0.286 0.480 0.740 1.059 1.410 1.746 2.018 2.182 2.219 2.134 1.955 1.719 1.461 1.208 0.977 0.776 0.609 0.472 0.363 0.279
[11,] 0.108 0.202 0.344 0.539 0.784 1.062 1.341 1.579 1.739 1.799 1.757 1.632 1.452 1.246 1.038 0.844 0.674 0.530 0.413 0.318 0.245
[12,] 0.076 0.142 0.245 0.387 0.570 0.783 1.002 1.197 1.337 1.401 1.385 1.299 1.166 1.008 0.844 0.690 0.553 0.437 0.340 0.263 0.203
[13,] 0.053 0.100 0.173 0.276 0.410 0.568 0.734 0.886 1.000 1.058 1.055 0.998 0.901 0.783 0.659 0.541 0.435 0.344 0.269 0.208 0.161
[14,] 0.037 0.070 0.122 0.196 0.292 0.408 0.531 0.645 0.733 0.781 0.784 0.746 0.677 0.591 0.499 0.411 0.331 0.262 0.205 0.159 0.123
[15,] 0.026 0.049 0.086 0.138 0.207 0.290 0.380 0.464 0.530 0.568 0.573 0.547 0.498 0.436 0.369 0.304 0.246 0.195 0.153 0.118 0.091
[16,] 0.018 0.034 0.060 0.097 0.146 0.205 0.270 0.331 0.380 0.408 0.413 0.395 0.361 0.317 0.269 0.222 0.179 0.142 0.111 0.086 0.067
[17,] 0.013 0.024 0.042 0.069 0.103 0.146 0.192 0.236 0.271 0.292 0.296 0.285 0.260 0.229 0.194 0.161 0.130 0.103 0.081 0.063 0.048

, , 8

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 2.662 2.020 1.522 1.146 0.862 0.648 0.487 0.366 0.275 0.207 0.155 0.116 0.087 0.066 0.049 0.037 0.028 0.021 0.016 0.012 0.009
 [2,] 1.882 1.874 1.754 1.563 1.338 1.112 0.901 0.718 0.563 0.437 0.337 0.258 0.196 0.149 0.113 0.085 0.064 0.048 0.036 0.027 0.021
 [3,] 1.321 1.608 1.789 1.843 1.779 1.629 1.429 1.210 0.996 0.803 0.636 0.498 0.385 0.296 0.226 0.172 0.131 0.099 0.075 0.056 0.043
 [4,] 0.926 1.302 1.656 1.926 2.071 2.082 1.975 1.785 1.550 1.303 1.067 0.857 0.677 0.528 0.408 0.313 0.239 0.182 0.138 0.104 0.079
 [5,] 0.648 1.011 1.423 1.825 2.154 2.359 2.416 2.333 2.142 1.886 1.603 1.325 1.072 0.851 0.667 0.518 0.398 0.305 0.232 0.176 0.134
 [6,] 0.454 0.761 1.156 1.600 2.040 2.407 2.645 2.723 2.644 2.441 2.159 1.843 1.529 1.240 0.987 0.775 0.602 0.464 0.355 0.271 0.207
 [7,] 0.317 0.561 0.899 1.321 1.790 2.249 2.627 2.866 2.934 2.838 2.612 2.305 1.965 1.628 1.319 1.049 0.823 0.639 0.492 0.377 0.289
 [8,] 0.222 0.407 0.679 1.041 1.480 1.954 2.402 2.754 2.955 2.984 2.852 2.599 2.275 1.926 1.588 1.281 1.016 0.795 0.616 0.474 0.364
 [9,] 0.155 0.292 0.501 0.793 1.168 1.603 2.052 2.451 2.738 2.869 2.836 2.661 2.388 2.064 1.729 1.414 1.133 0.894 0.697 0.539 0.416
[10,] 0.108 0.208 0.364 0.590 0.891 1.259 1.663 2.052 2.366 2.557 2.600 2.503 2.295 2.020 1.718 1.421 1.150 0.914 0.717 0.557 0.431
[11,] 0.076 0.147 0.261 0.430 0.663 0.957 1.294 1.636 1.935 2.143 2.230 2.191 2.046 1.828 1.575 1.316 1.074 0.859 0.677 0.528 0.410
[12,] 0.053 0.104 0.186 0.310 0.484 0.710 0.976 1.258 1.516 1.711 1.812 1.810 1.715 1.552 1.351 1.138 0.935 0.752 0.595 0.465 0.362
[13,] 0.037 0.073 0.132 0.221 0.349 0.518 0.721 0.941 1.150 1.317 1.414 1.430 1.370 1.251 1.098 0.932 0.769 0.621 0.494 0.387 0.301
[14,] 0.026 0.051 0.093 0.157 0.249 0.373 0.524 0.690 0.853 0.986 1.069 1.092 1.055 0.971 0.857 0.731 0.606 0.491 0.391 0.307 0.240
[15,] 0.018 0.036 0.065 0.111 0.177 0.266 0.376 0.499 0.621 0.724 0.791 0.813 0.790 0.731 0.648 0.555 0.461 0.375 0.299 0.235 0.184
[16,] 0.013 0.025 0.046 0.078 0.125 0.189 0.268 0.358 0.447 0.524 0.575 0.594 0.580 0.539 0.479 0.411 0.343 0.279 0.223 0.176 0.137
[17,] 0.009 0.018 0.032 0.055 0.089 0.134 0.191 0.256 0.321 0.377 0.416 0.431 0.422 0.393 0.351 0.301 0.252 0.205 0.164 0.129 0.101

, , 9

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 1.860 1.412 1.064 0.801 0.603 0.453 0.341 0.256 0.192 0.144 0.108 0.081 0.061 0.046 0.034 0.026 0.019 0.015 0.011 0.008 0.006
 [2,] 1.316 1.320 1.243 1.113 0.957 0.797 0.648 0.517 0.406 0.316 0.243 0.186 0.142 0.108 0.082 0.062 0.047 0.035 0.026 0.020 0.015
 [3,] 0.923 1.139 1.282 1.335 1.300 1.200 1.059 0.901 0.745 0.602 0.479 0.375 0.291 0.224 0.171 0.130 0.099 0.075 0.056 0.043 0.032
 [4,] 0.647 0.926 1.198 1.415 1.544 1.573 1.510 1.378 1.207 1.021 0.841 0.678 0.537 0.420 0.325 0.250 0.191 0.145 0.110 0.083 0.063
 [5,] 0.453 0.721 1.037 1.358 1.636 1.827 1.906 1.871 1.742 1.552 1.332 1.110 0.903 0.721 0.567 0.441 0.340 0.261 0.199 0.151 0.115
 [6,] 0.317 0.544 0.847 1.203 1.574 1.907 2.150 2.267 2.250 2.117 1.903 1.646 1.380 1.129 0.905 0.714 0.557 0.430 0.330 0.252 0.193
 [7,] 0.222 0.402 0.662 1.000 1.399 1.816 2.194 2.472 2.610 2.595 2.447 2.204 1.911 1.605 1.314 1.055 0.834 0.651 0.503 0.386 0.296
 [8,] 0.155 0.292 0.501 0.793 1.168 1.603 2.052 2.451 2.738 2.869 2.836 2.661 2.388 2.064 1.729 1.414 1.133 0.894 0.697 0.539 0.416
 [9,] 0.108 0.209 0.371 0.607 0.929 1.332 1.787 2.241 2.628 2.886 2.979 2.907 2.697 2.398 2.057 1.713 1.393 1.112 0.875 0.681 0.528
[10,] 0.076 0.149 0.270 0.453 0.713 1.056 1.469 1.915 2.338 2.672 2.866 2.895 2.770 2.528 2.217 1.880 1.551 1.253 0.994 0.779 0.607
[11,] 0.053 0.106 0.194 0.331 0.533 0.808 1.156 1.553 1.957 2.310 2.556 2.659 2.613 2.441 2.183 1.882 1.573 1.284 1.028 0.810 0.635
[12,] 0.037 0.074 0.138 0.239 0.390 0.603 0.880 1.209 1.561 1.889 2.144 2.284 2.294 2.186 1.987 1.737 1.469 1.210 0.976 0.773 0.608
[13,] 0.026 0.052 0.098 0.171 0.282 0.441 0.653 0.913 1.200 1.481 1.713 1.859 1.901 1.839 1.695 1.498 1.279 1.061 0.861 0.686 0.541
[14,] 0.018 0.037 0.069 0.121 0.202 0.319 0.477 0.674 0.898 1.124 1.319 1.452 1.504 1.473 1.372 1.224 1.053 0.879 0.716 0.573 0.453
[15,] 0.013 0.026 0.049 0.086 0.143 0.228 0.343 0.490 0.659 0.833 0.988 1.099 1.150 1.137 1.067 0.958 0.829 0.695 0.569 0.456 0.362
[16,] 0.009 0.018 0.034 0.060 0.101 0.162 0.245 0.352 0.477 0.607 0.726 0.814 0.857 0.853 0.805 0.727 0.632 0.532 0.436 0.351 0.279
[17,] 0.006 0.013 0.024 0.043 0.072 0.115 0.175 0.252 0.343 0.439 0.528 0.595 0.630 0.630 0.598 0.541 0.472 0.398 0.327 0.263 0.210

, , 10

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 1.300 0.986 0.744 0.560 0.421 0.317 0.238 0.179 0.134 0.101 0.076 0.057 0.043 0.032 0.024 0.018 0.014 0.010 0.008 0.006 0.004
 [2,] 0.919 0.927 0.877 0.788 0.680 0.568 0.462 0.369 0.290 0.226 0.174 0.133 0.102 0.077 0.059 0.044 0.033 0.025 0.019 0.014 0.011
 [3,] 0.645 0.803 0.912 0.956 0.938 0.870 0.771 0.659 0.546 0.443 0.352 0.276 0.215 0.165 0.126 0.096 0.073 0.055 0.042 0.032 0.024
 [4,] 0.452 0.655 0.858 1.025 1.130 1.162 1.125 1.034 0.911 0.774 0.640 0.518 0.411 0.322 0.250 0.192 0.147 0.112 0.085 0.064 0.049
 [5,] 0.317 0.511 0.746 0.992 1.213 1.375 1.453 1.444 1.359 1.221 1.056 0.885 0.724 0.580 0.458 0.357 0.276 0.211 0.161 0.123 0.093
 [6,] 0.222 0.387 0.612 0.885 1.180 1.459 1.677 1.801 1.817 1.735 1.579 1.381 1.168 0.962 0.775 0.614 0.481 0.372 0.286 0.219 0.167
 [7,] 0.155 0.286 0.480 0.740 1.059 1.410 1.746 2.018 2.182 2.219 2.134 1.955 1.719 1.461 1.208 0.977 0.776 0.609 0.472 0.363 0.279
 [8,] 0.108 0.208 0.364 0.590 0.891 1.259 1.663 2.052 2.366 2.557 2.600 2.503 2.295 2.020 1.718 1.421 1.150 0.914 0.717 0.557 0.431
 [9,] 0.076 0.149 0.270 0.453 0.713 1.056 1.469 1.915 2.338 2.672 2.866 2.895 2.770 2.528 2.217 1.880 1.551 1.253 0.994 0.779 0.607
[10,] 0.053 0.106 0.196 0.338 0.550 0.843 1.222 1.666 2.133 2.558 2.877 3.039 3.029 2.865 2.589 2.251 1.896 1.556 1.252 0.990 0.778
[11,] 0.037 0.075 0.141 0.248 0.412 0.649 0.970 1.369 1.821 2.273 2.661 2.923 3.022 2.953 2.746 2.446 2.101 1.752 1.427 1.141 0.903
[12,] 0.026 0.053 0.101 0.179 0.302 0.486 0.743 1.077 1.475 1.900 2.299 2.608 2.781 2.796 2.667 2.427 2.123 1.797 1.481 1.195 0.952
[13,] 0.018 0.037 0.071 0.128 0.219 0.357 0.555 0.820 1.148 1.515 1.880 2.189 2.393 2.464 2.399 2.224 1.976 1.695 1.411 1.148 0.921
[14,] 0.013 0.026 0.050 0.091 0.157 0.258 0.406 0.609 0.866 1.165 1.473 1.750 1.951 2.046 2.027 1.908 1.717 1.488 1.250 1.024 0.826
[15,] 0.009 0.018 0.035 0.064 0.111 0.185 0.293 0.444 0.640 0.871 1.118 1.348 1.526 1.623 1.629 1.552 1.411 1.233 1.043 0.859 0.696
[16,] 0.006 0.013 0.025 0.045 0.079 0.131 0.210 0.320 0.465 0.639 0.829 1.010 1.156 1.243 1.261 1.211 1.109 0.976 0.830 0.687 0.558
[17,] 0.004 0.009 0.018 0.032 0.056 0.093 0.150 0.230 0.336 0.465 0.607 0.746 0.860 0.932 0.952 0.921 0.848 0.750 0.640 0.531 0.433

, , 11

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 0.908 0.689 0.520 0.391 0.294 0.221 0.166 0.125 0.094 0.071 0.053 0.040 0.030 0.022 0.017 0.013 0.009 0.007 0.005 0.004 0.003
 [2,] 0.642 0.650 0.617 0.556 0.480 0.402 0.328 0.262 0.206 0.161 0.124 0.095 0.072 0.055 0.042 0.031 0.024 0.018 0.013 0.010 0.008
 [3,] 0.451 0.565 0.645 0.680 0.670 0.624 0.555 0.475 0.395 0.320 0.255 0.201 0.156 0.120 0.092 0.070 0.053 0.040 0.030 0.023 0.017
 [4,] 0.316 0.462 0.610 0.734 0.816 0.845 0.823 0.760 0.672 0.574 0.476 0.385 0.307 0.241 0.187 0.144 0.110 0.084 0.064 0.048 0.037
 [5,] 0.221 0.361 0.532 0.715 0.884 1.012 1.081 1.083 1.027 0.930 0.808 0.681 0.558 0.449 0.355 0.277 0.214 0.165 0.126 0.096 0.073
 [6,] 0.155 0.273 0.438 0.641 0.867 1.087 1.268 1.380 1.410 1.362 1.252 1.103 0.939 0.778 0.629 0.501 0.393 0.305 0.235 0.180 0.138
 [7,] 0.108 0.202 0.344 0.539 0.784 1.062 1.341 1.579 1.739 1.799 1.757 1.632 1.452 1.246 1.038 0.844 0.674 0.530 0.413 0.318 0.245
 [8,] 0.076 0.147 0.261 0.430 0.663 0.957 1.294 1.636 1.935 2.143 2.230 2.191 2.046 1.828 1.575 1.316 1.074 0.859 0.677 0.528 0.410
 [9,] 0.053 0.106 0.194 0.331 0.533 0.808 1.156 1.553 1.957 2.310 2.556 2.659 2.613 2.441 2.183 1.882 1.573 1.284 1.028 0.810 0.635
[10,] 0.037 0.075 0.141 0.248 0.412 0.649 0.970 1.369 1.821 2.273 2.661 2.923 3.022 2.953 2.746 2.446 2.101 1.752 1.427 1.141 0.903
[11,] 0.026 0.053 0.102 0.182 0.309 0.502 0.775 1.138 1.580 2.067 2.541 2.931 3.174 3.239 3.130 2.881 2.545 2.172 1.801 1.460 1.169
[12,] 0.018 0.038 0.073 0.132 0.228 0.377 0.597 0.902 1.296 1.761 2.253 2.709 3.056 3.242 3.246 3.083 2.798 2.441 2.062 1.696 1.373
[13,] 0.013 0.026 0.051 0.094 0.165 0.277 0.447 0.691 1.018 1.424 1.882 2.339 2.731 2.994 3.092 3.020 2.807 2.500 2.148 1.791 1.465
[14,] 0.009 0.019 0.036 0.067 0.118 0.201 0.328 0.516 0.774 1.107 1.499 1.912 2.294 2.584 2.737 2.737 2.598 2.355 2.054 1.734 1.432
[15,] 0.006 0.013 0.026 0.047 0.084 0.144 0.238 0.377 0.575 0.835 1.151 1.499 1.835 2.111 2.283 2.326 2.246 2.067 1.825 1.556 1.296
[16,] 0.004 0.009 0.018 0.033 0.060 0.102 0.170 0.273 0.419 0.616 0.861 1.137 1.415 1.654 1.817 1.879 1.838 1.711 1.526 1.313 1.100
[17,] 0.003 0.006 0.013 0.024 0.042 0.073 0.122 0.196 0.304 0.450 0.635 0.847 1.066 1.261 1.401 1.465 1.448 1.360 1.222 1.058 0.891

, , 12

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 0.634 0.482 0.363 0.273 0.206 0.155 0.116 0.087 0.066 0.049 0.037 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004 0.003 0.002
 [2,] 0.449 0.456 0.433 0.391 0.338 0.283 0.231 0.185 0.146 0.113 0.088 0.067 0.051 0.039 0.029 0.022 0.017 0.013 0.010 0.007 0.005
 [3,] 0.315 0.397 0.455 0.481 0.475 0.444 0.396 0.339 0.282 0.230 0.183 0.144 0.112 0.086 0.066 0.050 0.038 0.029 0.022 0.016 0.012
 [4,] 0.221 0.325 0.431 0.522 0.583 0.607 0.594 0.550 0.488 0.418 0.347 0.282 0.224 0.176 0.137 0.106 0.081 0.062 0.047 0.035 0.027
 [5,] 0.155 0.254 0.377 0.511 0.636 0.734 0.789 0.796 0.759 0.690 0.603 0.509 0.419 0.337 0.267 0.209 0.162 0.124 0.095 0.072 0.055
 [6,] 0.108 0.192 0.311 0.460 0.628 0.795 0.937 1.030 1.063 1.035 0.958 0.849 0.727 0.604 0.491 0.391 0.307 0.239 0.184 0.141 0.108
 [7,] 0.076 0.142 0.245 0.387 0.570 0.783 1.002 1.197 1.337 1.401 1.385 1.299 1.166 1.008 0.844 0.690 0.553 0.437 0.340 0.263 0.203
 [8,] 0.053 0.104 0.186 0.310 0.484 0.710 0.976 1.258 1.516 1.711 1.812 1.810 1.715 1.552 1.351 1.138 0.935 0.752 0.595 0.465 0.362
 [9,] 0.037 0.074 0.138 0.239 0.390 0.603 0.880 1.209 1.561 1.889 2.144 2.284 2.294 2.186 1.987 1.737 1.469 1.210 0.976 0.773 0.608
[10,] 0.026 0.053 0.101 0.179 0.302 0.486 0.743 1.077 1.475 1.900 2.299 2.608 2.781 2.796 2.667 2.427 2.123 1.797 1.481 1.195 0.952
[11,] 0.018 0.038 0.073 0.132 0.228 0.377 0.597 0.902 1.296 1.761 2.253 2.709 3.056 3.242 3.246 3.083 2.798 2.441 2.062 1.696 1.373
[12,] 0.013 0.027 0.052 0.095 0.168 0.284 0.461 0.720 1.074 1.523 2.043 2.582 3.065 3.417 3.585 3.553 3.345 3.012 2.611 2.194 1.805
[13,] 0.009 0.019 0.037 0.068 0.122 0.209 0.347 0.554 0.851 1.247 1.737 2.287 2.835 3.302 3.613 3.724 3.633 3.373 3.001 2.576 2.155
[14,] 0.006 0.013 0.026 0.048 0.087 0.152 0.255 0.415 0.651 0.978 1.403 1.908 2.449 2.958 3.356 3.581 3.606 3.444 3.139 2.750 2.339
[15,] 0.004 0.009 0.018 0.034 0.062 0.109 0.185 0.305 0.485 0.743 1.089 1.519 2.004 2.491 2.911 3.197 3.307 3.236 3.013 2.688 2.320
[16,] 0.003 0.006 0.013 0.024 0.044 0.077 0.133 0.220 0.355 0.551 0.821 1.166 1.571 1.997 2.389 2.684 2.838 2.834 2.686 2.433 2.126
[17,] 0.002 0.005 0.009 0.017 0.031 0.055 0.095 0.159 0.257 0.404 0.608 0.876 1.198 1.549 1.885 2.155 2.318 2.351 2.260 2.072 1.830

, , 13

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 0.443 0.337 0.254 0.191 0.144 0.108 0.081 0.061 0.046 0.034 0.026 0.019 0.015 0.011 0.008 0.006 0.005 0.003 0.003 0.002 0.001
 [2,] 0.314 0.319 0.304 0.274 0.238 0.199 0.163 0.130 0.103 0.080 0.062 0.047 0.036 0.027 0.021 0.016 0.012 0.009 0.007 0.005 0.004
 [3,] 0.220 0.278 0.320 0.339 0.336 0.314 0.280 0.241 0.201 0.163 0.130 0.102 0.080 0.061 0.047 0.036 0.027 0.021 0.016 0.012 0.009
 [4,] 0.154 0.228 0.304 0.369 0.414 0.432 0.424 0.394 0.351 0.301 0.250 0.203 0.162 0.127 0.099 0.076 0.058 0.045 0.034 0.026 0.019
 [5,] 0.108 0.178 0.266 0.363 0.454 0.527 0.569 0.577 0.552 0.504 0.441 0.374 0.308 0.248 0.197 0.154 0.119 0.092 0.070 0.053 0.041
 [6,] 0.076 0.135 0.220 0.327 0.450 0.574 0.681 0.754 0.783 0.767 0.714 0.636 0.546 0.456 0.371 0.296 0.233 0.182 0.140 0.107 0.082
 [7,] 0.053 0.100 0.173 0.276 0.410 0.568 0.734 0.886 1.000 1.058 1.055 0.998 0.901 0.783 0.659 0.541 0.435 0.344 0.269 0.208 0.161
 [8,] 0.037 0.073 0.132 0.221 0.349 0.518 0.721 0.941 1.150 1.317 1.414 1.430 1.370 1.251 1.098 0.932 0.769 0.621 0.494 0.387 0.301
 [9,] 0.026 0.052 0.098 0.171 0.282 0.441 0.653 0.913 1.200 1.481 1.713 1.859 1.901 1.839 1.695 1.498 1.279 1.061 0.861 0.686 0.541
[10,] 0.018 0.037 0.071 0.128 0.219 0.357 0.555 0.820 1.148 1.515 1.880 2.189 2.393 2.464 2.399 2.224 1.976 1.695 1.411 1.148 0.921
[11,] 0.013 0.026 0.051 0.094 0.165 0.277 0.447 0.691 1.018 1.424 1.882 2.339 2.731 2.994 3.092 3.020 2.807 2.500 2.148 1.791 1.465
[12,] 0.009 0.019 0.037 0.068 0.122 0.209 0.347 0.554 0.851 1.247 1.737 2.287 2.835 3.302 3.613 3.724 3.633 3.373 3.001 2.576 2.155
[13,] 0.006 0.013 0.026 0.049 0.088 0.154 0.261 0.428 0.677 1.031 1.498 2.070 2.702 3.323 3.839 4.169 4.267 4.135 3.819 3.382 2.902
[14,] 0.004 0.009 0.018 0.035 0.063 0.112 0.192 0.321 0.521 0.814 1.224 1.757 2.394 3.083 3.736 4.254 4.556 4.604 4.411 4.033 3.552
[15,] 0.003 0.006 0.013 0.025 0.045 0.080 0.140 0.236 0.389 0.622 0.958 1.417 1.998 2.671 3.369 3.996 4.456 4.677 4.638 4.370 3.948
[16,] 0.002 0.005 0.009 0.017 0.032 0.057 0.100 0.171 0.286 0.463 0.727 1.099 1.591 2.190 2.852 3.498 4.035 4.375 4.473 4.332 4.005
[17,] 0.001 0.003 0.006 0.012 0.023 0.041 0.072 0.123 0.208 0.340 0.541 0.832 1.227 1.727 2.304 2.902 3.438 3.828 4.013 3.976 3.750

, , 14

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18] [,19] [,20] [,21]
 [1,] 0.310 0.235 0.177 0.134 0.100 0.076 0.057 0.043 0.032 0.024 0.018 0.014 0.010 0.008 0.006 0.004 0.003 0.002 0.002 0.001 0.001
 [2,] 0.219 0.223 0.213 0.192 0.167 0.140 0.114 0.091 0.072 0.056 0.043 0.033 0.025 0.019 0.015 0.011 0.008 0.006 0.005 0.004 0.003
 [3,] 0.154 0.195 0.224 0.238 0.236 0.222 0.198 0.170 0.142 0.115 0.092 0.073 0.056 0.044 0.033 0.025 0.019 0.015 0.011 0.008 0.006
 [4,] 0.108 0.160 0.214 0.260 0.293 0.306 0.301 0.281 0.250 0.214 0.179 0.145 0.116 0.091 0.071 0.055 0.042 0.032 0.024 0.018 0.014
 [5,] 0.076 0.125 0.187 0.256 0.322 0.375 0.406 0.413 0.397 0.363 0.319 0.270 0.223 0.180 0.143 0.112 0.087 0.067 0.051 0.039 0.030
 [6,] 0.053 0.095 0.155 0.231 0.320 0.410 0.490 0.545 0.569 0.559 0.523 0.467 0.402 0.336 0.274 0.219 0.173 0.135 0.104 0.080 0.061
 [7,] 0.037 0.070 0.122 0.196 0.292 0.408 0.531 0.645 0.733 0.781 0.784 0.746 0.677 0.591 0.499 0.411 0.331 0.262 0.205 0.159 0.123
 [8,] 0.026 0.051 0.093 0.157 0.249 0.373 0.524 0.690 0.853 0.986 1.069 1.092 1.055 0.971 0.857 0.731 0.606 0.491 0.391 0.307 0.240
 [9,] 0.018 0.037 0.069 0.121 0.202 0.319 0.477 0.674 0.898 1.124 1.319 1.452 1.504 1.473 1.372 1.224 1.053 0.879 0.716 0.573 0.453
[10,] 0.013 0.026 0.050 0.091 0.157 0.258 0.406 0.609 0.866 1.165 1.473 1.750 1.951 2.046 2.027 1.908 1.717 1.488 1.250 1.024 0.826
[11,] 0.009 0.019 0.036 0.067 0.118 0.201 0.328 0.516 0.774 1.107 1.499 1.912 2.294 2.584 2.737 2.737 2.598 2.355 2.054 1.734 1.432
[12,] 0.006 0.013 0.026 0.048 0.087 0.152 0.255 0.415 0.651 0.978 1.403 1.908 2.449 2.958 3.356 3.581 3.606 3.444 3.139 2.750 2.339
[13,] 0.004 0.009 0.018 0.035 0.063 0.112 0.192 0.321 0.521 0.814 1.224 1.757 2.394 3.083 3.736 4.254 4.556 4.604 4.411 4.033 3.552
[14,] 0.003 0.006 0.013 0.025 0.045 0.081 0.142 0.242 0.401 0.647 1.009 1.512 2.167 2.949 3.793 4.593 5.227 5.595 5.651 5.413 4.958
[15,] 0.002 0.005 0.009 0.017 0.032 0.058 0.103 0.178 0.301 0.496 0.796 1.234 1.839 2.621 3.547 4.535 5.459 6.176 6.575 6.607 6.309
[16,] 0.001 0.003 0.006 0.012 0.023 0.042 0.074 0.129 0.221 0.371 0.607 0.965 1.484 2.193 3.093 4.141 5.234 6.225 6.960 7.326 7.290
[17,] 0.001 0.002 0.005 0.009 0.016 0.030 0.053 0.093 0.161 0.273 0.453 0.735 1.156 1.756 2.558 3.552 4.672 5.795 6.762 7.419 7.669

, , 15

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17] [,18]  [,19]  [,20]  [,21]
 [1,] 0.216 0.164 0.124 0.093 0.070 0.053 0.040 0.030 0.022 0.017 0.013 0.009 0.007 0.005 0.004 0.003 0.002 0.002  0.001  0.001  0.001
 [2,] 0.153 0.156 0.149 0.135 0.117 0.098 0.080 0.064 0.051 0.039 0.030 0.023 0.018 0.014 0.010 0.008 0.006 0.004  0.003  0.002  0.002
 [3,] 0.108 0.136 0.157 0.167 0.166 0.156 0.139 0.120 0.100 0.081 0.065 0.051 0.040 0.031 0.024 0.018 0.014 0.010  0.008  0.006  0.004
 [4,] 0.075 0.112 0.150 0.183 0.206 0.216 0.213 0.198 0.177 0.152 0.127 0.103 0.082 0.065 0.050 0.039 0.030 0.023  0.017  0.013  0.010
 [5,] 0.053 0.088 0.132 0.180 0.227 0.265 0.288 0.294 0.283 0.259 0.228 0.193 0.160 0.129 0.103 0.080 0.062 0.048  0.037  0.028  0.021
 [6,] 0.037 0.066 0.109 0.163 0.226 0.291 0.349 0.390 0.408 0.403 0.377 0.338 0.292 0.244 0.199 0.160 0.126 0.098  0.076  0.058  0.045
 [7,] 0.026 0.049 0.086 0.138 0.207 0.290 0.380 0.464 0.530 0.568 0.573 0.547 0.498 0.436 0.369 0.304 0.246 0.195  0.153  0.118  0.091
 [8,] 0.018 0.036 0.065 0.111 0.177 0.266 0.376 0.499 0.621 0.724 0.791 0.813 0.790 0.731 0.648 0.555 0.461 0.375  0.299  0.235  0.184
 [9,] 0.013 0.026 0.049 0.086 0.143 0.228 0.343 0.490 0.659 0.833 0.988 1.099 1.150 1.137 1.067 0.958 0.829 0.695  0.569  0.456  0.362
[10,] 0.009 0.018 0.035 0.064 0.111 0.185 0.293 0.444 0.640 0.871 1.118 1.348 1.526 1.623 1.629 1.552 1.411 1.233  1.043  0.859  0.696
[11,] 0.006 0.013 0.026 0.047 0.084 0.144 0.238 0.377 0.575 0.835 1.151 1.499 1.835 2.111 2.283 2.326 2.246 2.067  1.825  1.556  1.296
[12,] 0.004 0.009 0.018 0.034 0.062 0.109 0.185 0.305 0.485 0.743 1.089 1.519 2.004 2.491 2.911 3.197 3.307 3.236  3.013  2.688  2.320
[13,] 0.003 0.006 0.013 0.025 0.045 0.080 0.140 0.236 0.389 0.622 0.958 1.417 1.998 2.671 3.369 3.996 4.456 4.677  4.638  4.370  3.948
[14,] 0.002 0.005 0.009 0.017 0.032 0.058 0.103 0.178 0.301 0.496 0.796 1.234 1.839 2.621 3.547 4.535 5.459 6.176  6.575  6.607  6.309
[15,] 0.001 0.003 0.006 0.012 0.023 0.042 0.075 0.131 0.226 0.382 0.631 1.015 1.583 2.380 3.425 4.689 6.073 7.412  8.505  9.179  9.346
[16,] 0.001 0.002 0.005 0.009 0.016 0.030 0.054 0.095 0.166 0.286 0.483 0.799 1.292 2.026 3.068 4.456 6.166 8.084  9.995 11.622 12.694
[17,] 0.001 0.002 0.003 0.006 0.012 0.021 0.039 0.069 0.121 0.211 0.362 0.612 1.015 1.644 2.592 3.948 5.777 8.065 10.680 13.346 15.651

, , 16

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17]  [,18]  [,19]  [,20]  [,21]
 [1,] 0.151 0.115 0.087 0.065 0.049 0.037 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004 0.003 0.002 0.002  0.001  0.001  0.001  0.001
 [2,] 0.107 0.109 0.104 0.094 0.082 0.069 0.056 0.045 0.035 0.028 0.021 0.016 0.012 0.009 0.007 0.005 0.004  0.003  0.002  0.002  0.001
 [3,] 0.075 0.095 0.110 0.117 0.117 0.109 0.098 0.084 0.070 0.057 0.046 0.036 0.028 0.022 0.017 0.013 0.010  0.007  0.005  0.004  0.003
 [4,] 0.053 0.078 0.105 0.128 0.145 0.152 0.150 0.140 0.125 0.107 0.089 0.073 0.058 0.046 0.036 0.027 0.021  0.016  0.012  0.009  0.007
 [5,] 0.037 0.061 0.092 0.127 0.160 0.187 0.204 0.208 0.200 0.184 0.162 0.137 0.114 0.092 0.073 0.057 0.044  0.034  0.026  0.020  0.015
 [6,] 0.026 0.047 0.076 0.115 0.160 0.206 0.247 0.277 0.291 0.287 0.270 0.242 0.209 0.175 0.143 0.115 0.091  0.071  0.055  0.042  0.032
 [7,] 0.018 0.034 0.060 0.097 0.146 0.205 0.270 0.331 0.380 0.408 0.413 0.395 0.361 0.317 0.269 0.222 0.179  0.142  0.111  0.086  0.067
 [8,] 0.013 0.025 0.046 0.078 0.125 0.189 0.268 0.358 0.447 0.524 0.575 0.594 0.580 0.539 0.479 0.411 0.343  0.279  0.223  0.176  0.137
 [9,] 0.009 0.018 0.034 0.060 0.101 0.162 0.245 0.352 0.477 0.607 0.726 0.814 0.857 0.853 0.805 0.727 0.632  0.532  0.436  0.351  0.279
[10,] 0.006 0.013 0.025 0.045 0.079 0.131 0.210 0.320 0.465 0.639 0.829 1.010 1.156 1.243 1.261 1.211 1.109  0.976  0.830  0.687  0.558
[11,] 0.004 0.009 0.018 0.033 0.060 0.102 0.170 0.273 0.419 0.616 0.861 1.137 1.415 1.654 1.817 1.879 1.838  1.711  1.526  1.313  1.100
[12,] 0.003 0.006 0.013 0.024 0.044 0.077 0.133 0.220 0.355 0.551 0.821 1.166 1.571 1.997 2.389 2.684 2.838  2.834  2.686  2.433  2.126
[13,] 0.002 0.005 0.009 0.017 0.032 0.057 0.100 0.171 0.286 0.463 0.727 1.099 1.591 2.190 2.852 3.498 4.035  4.375  4.473  4.332  4.005
[14,] 0.001 0.003 0.006 0.012 0.023 0.042 0.074 0.129 0.221 0.371 0.607 0.965 1.484 2.193 3.093 4.141 5.234  6.225  6.960  7.326  7.290
[15,] 0.001 0.002 0.005 0.009 0.016 0.030 0.054 0.095 0.166 0.286 0.483 0.799 1.292 2.026 3.068 4.456 6.166  8.084  9.995 11.622 12.694
[16,] 0.001 0.002 0.003 0.006 0.012 0.021 0.039 0.069 0.122 0.214 0.371 0.633 1.063 1.750 2.812 4.388 6.607  9.538 13.115 17.075 20.887
[17,] 0.001 0.001 0.002 0.004 0.008 0.015 0.028 0.050 0.089 0.158 0.279 0.486 0.840 1.436 2.419 4.005 6.488 10.230 15.609 22.892 31.860

, , 17

       [,1]  [,2]  [,3]  [,4]  [,5]  [,6]  [,7]  [,8]  [,9] [,10] [,11] [,12] [,13] [,14] [,15] [,16] [,17]  [,18]  [,19]  [,20]  [,21]
 [1,] 0.106 0.081 0.061 0.046 0.034 0.026 0.019 0.015 0.011 0.008 0.006 0.005 0.003 0.003 0.002 0.001 0.001  0.001  0.001  0.000  0.000
 [2,] 0.075 0.077 0.073 0.066 0.058 0.048 0.039 0.032 0.025 0.019 0.015 0.012 0.009 0.007 0.005 0.004 0.003  0.002  0.002  0.001  0.001
 [3,] 0.053 0.067 0.077 0.083 0.082 0.077 0.069 0.059 0.050 0.040 0.032 0.025 0.020 0.015 0.012 0.009 0.007  0.005  0.004  0.003  0.002
 [4,] 0.037 0.055 0.074 0.091 0.102 0.107 0.106 0.099 0.088 0.076 0.063 0.052 0.041 0.032 0.025 0.019 0.015  0.011  0.009  0.007  0.005
 [5,] 0.026 0.043 0.065 0.089 0.113 0.132 0.144 0.147 0.142 0.131 0.115 0.098 0.081 0.065 0.052 0.041 0.032  0.024  0.019  0.014  0.011
 [6,] 0.018 0.033 0.054 0.081 0.113 0.146 0.175 0.197 0.207 0.205 0.193 0.173 0.150 0.125 0.103 0.082 0.065  0.051  0.039  0.030  0.023
 [7,] 0.013 0.024 0.042 0.069 0.103 0.146 0.192 0.236 0.271 0.292 0.296 0.285 0.260 0.229 0.194 0.161 0.130  0.103  0.081  0.063  0.048
 [8,] 0.009 0.018 0.032 0.055 0.089 0.134 0.191 0.256 0.321 0.377 0.416 0.431 0.422 0.393 0.351 0.301 0.252  0.205  0.164  0.129  0.101
 [9,] 0.006 0.013 0.024 0.043 0.072 0.115 0.175 0.252 0.343 0.439 0.528 0.595 0.630 0.630 0.598 0.541 0.472  0.398  0.327  0.263  0.210
[10,] 0.004 0.009 0.018 0.032 0.056 0.093 0.150 0.230 0.336 0.465 0.607 0.746 0.860 0.932 0.952 0.921 0.848  0.750  0.640  0.531  0.433
[11,] 0.003 0.006 0.013 0.024 0.042 0.073 0.122 0.196 0.304 0.450 0.635 0.847 1.066 1.261 1.401 1.465 1.448  1.360  1.222  1.058  0.891
[12,] 0.002 0.005 0.009 0.017 0.031 0.055 0.095 0.159 0.257 0.404 0.608 0.876 1.198 1.549 1.885 2.155 2.318  2.351  2.260  2.072  1.830
[13,] 0.001 0.003 0.006 0.012 0.023 0.041 0.072 0.123 0.208 0.340 0.541 0.832 1.227 1.727 2.304 2.902 3.438  3.828  4.013  3.976  3.750
[14,] 0.001 0.002 0.005 0.009 0.016 0.030 0.053 0.093 0.161 0.273 0.453 0.735 1.156 1.756 2.558 3.552 4.672  5.795  6.762  7.419  7.669
[15,] 0.001 0.002 0.003 0.006 0.012 0.021 0.039 0.069 0.121 0.211 0.362 0.612 1.015 1.644 2.592 3.948 5.777  8.065 10.680 13.346 15.651
[16,] 0.001 0.001 0.002 0.004 0.008 0.015 0.028 0.050 0.089 0.158 0.279 0.486 0.840 1.436 2.419 4.005 6.488 10.230 15.609 22.892 31.860
[17,] 0.000 0.001 0.002 0.003 0.006 0.011 0.020 0.036 0.065 0.117 0.210 0.375 0.668 1.188 2.112 3.750 6.647 11.769 20.806 36.713 63.909

> 
> ## dim = 4 --- fine as long we are "out of corners"
> um4 <- Exp.grid(u1=u1, u2=u2, u3=u1, u4=rev(u2))
> dcop <- dCopula(um4, frankCopula(param=theta.fr, dim = 4))
> stopifnot(dcop >= 0,# no NA here
+ 	  all.equal(dcop, copFrank@dacopula(um4, theta = theta.fr)))
> ## round(array(dcop, c(n1,n2,n1,n2)), 3)
> 
> showProc.time()
Time (user system elapsed): 1.52 0.03 1.55 
> 
> 
> ### --- now look at dCopula() and pCopula() for *all* copulas:
> 
> ##' u "mostly" in [0,1] .. with exceptions that are even NA, NaN
> mku <- function(n, fr = 1/20, sd = 1/10) {
+     r <- runif(n) + rnorm(n, sd=sd)
+     r[sample(n, ceiling(fr*n))] <- c(NA,NaN)
+     r
+ }
> 
> ## This is from ./moments.R --- keep in sync! ---
> tau.s <- c(       -.1, 0, 0.05805, (1:2)/9, 0.3)
> names(tau.s) <- paste0("tau=", sub("0[.]", ".", formatC(tau.s)))
> suppressWarnings( # warnings about ... "tau must be >= 0.
+     tTau <- sapply(tau.s, function(tau) vapply(copObs, iTau, numeric(1), tau = tau))
+ ) # tTau is printed in moments.R
> 
> suppressWarnings(RNGversion("3.5.0")) ; set.seed(12)
> u <- matrix(mku(1000), ncol= 2)# d = 2 required in the copula objects
> 
> ##' == function(u) copula:::outside.01(u, strictly = FALSE) --- used in dCopula()
> u.outside.01 <- function(u) apply(u, 1, function(x) any(x <= 0, 1 <= x))
> u.out <- u.outside.01(u) # has NAs
> u.ina <- apply(u, 1, anyNA)
> u.iNA <- is.na(u.out)
> (tuna <- table(u.ina, u.out, exclude={}))
       u.out
u.ina   FALSE TRUE <NA>
  FALSE   384   67    0
  TRUE      0    7   42
> ##        u.out
> ## u.ina   FALSE TRUE <NA>
> ##   FALSE   384   67    0
> ##   TRUE      0    7   42
> stopifnot(identical(c(384L, 0L, 67L, 7L, 0L, 42L), c(tuna)))
> u.OUT <- u.out & !u.ina  # no  NAs
> ## NB:  f(u) = 0 if some u_j is outside (0,1) - even when another u_k is NA
> 
> copsT <- lapply(setNames(,names(copObs)), # so have *named* list
+                  function(cN)
+                      lapply(unique(tTau[cN,]), function(th) setPar(copObs[[cN]], th)))
> copsTh <- unlist(copsT, recursive=FALSE)
> ## --> a list of 72 parametrized copulas
> ## ...  currently one of them not being "valid":
> okC <- lapply(copsTh, validObject, test=TRUE)
> (inOk <- which(notOk <- !sapply(okC, isTRUE)))
joeCopula2 
        43 
> ## joeCopula2
> ##         43
> if(length(inOk)) { # maybe no longer in future
+     stopifnot(identical(inOk, c(joeCopula2 = 43L)))
+     jp <- copsTh[[inOk]]@parameters
+     print(jp - 1) # -3.41e-10
+     stopifnot(all.equal(jp, 1))
+     ## now fix it up:
+     copsTh[[inOk]]@parameters <- 1
+     validObject(copsTh[[inOk]])
+ }
[1] -3.411909e-10
[1] TRUE
> ## Now, *all* are validObject :
> stopifnot( vapply(copsTh, validObject, logical(1), test=TRUE) )
> 
> 
> ## now using  ( u[], u.ina, u.OUT ):
> op <- options(warn = 1) # immediate
> ## ==>  *all* Huessler-Reiss  produce warnings  such as  << FIXME eventually
> ##  Warning in log(u1p/u2p) : NaNs produced
> ## all other copulas give no warning !
> errH <- if(getRversion() < "3.5.0") identity else function(e) {e$call <- NULL ;  e}
> copChk <- lapply(copsTh, function(cop) {
+     show(cop)
+     fu <- dCopula(u, cop)
+     Fu <- pCopula(u, cop)
+     lfu <- dCopula(u, cop, log=TRUE)
+     ## lfu. <- lfu[!u.ina]
+     cat("-----------\n")
+     tryCatch(
+         stopifnot(fu[u.OUT] == 0,
+ 		  all.equal(fu, exp(lfu), tolerance = 1e-15),
+                   0 <= fu[!u.ina],
+                   0 <= Fu[!u.ina], Fu[!u.ina] <= 1,
+                   is.finite(lfu[!u.ina]) | lfu[!u.ina] == -Inf,
+                   is.na(fu[u.iNA]), is.na(Fu[u.ina]), is.na(lfu[u.iNA]))
+        , error = errH)
+     })
Farlie-Gumbel-Morgenstern (FGM) copula, dim. d = 2
 param.: -0.45 
Dimension:  2 
Parameters:
  alpha1.2   = -0.45
-----------
Farlie-Gumbel-Morgenstern (FGM) copula, dim. d = 2
 param.: 0 
Dimension:  2 
Parameters:
  alpha1.2   = 0
-----------
Farlie-Gumbel-Morgenstern (FGM) copula, dim. d = 2
 param.: 0.261 
Dimension:  2 
Parameters:
  alpha1.2   = 0.261225
-----------
Farlie-Gumbel-Morgenstern (FGM) copula, dim. d = 2
 param.: 0.5 
Dimension:  2 
Parameters:
  alpha1.2   = 0.5
-----------
Farlie-Gumbel-Morgenstern (FGM) copula, dim. d = 2
 param.: 1 
Dimension:  2 
Parameters:
  alpha1.2   = 1
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.6395724
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.001153
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.297035
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.643074
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 2.754562
-----------
plackettCopula copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 3.986747
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = -0.1564345
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.09105842
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.1736482
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.3420201
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
Normal copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.4539905
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
Warning in probval.Miwa(algorithm, n, df, lower, upper, infin, corr, delta) :
  Approximating +/-Inf by +/-1000
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = -0.1564345
  df      =  4.0000000
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0
  df      = 4
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.09105842
  df      = 4.00000000
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.1736482
  df      = 4.0000000
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.3420201
  df      = 4.0000000
-----------
t-copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.4539905
  df      = 4.0000000
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = -0.1818182
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.1232549
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.25
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.5714286
-----------
Clayton copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.8571429
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = -0.9073675
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.5238811
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.010132
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 2.084387
-----------
Frank copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 2.917434
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = -0.5030297
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.368215e-13
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.2445963
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.4404231
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.7715254
-----------
AMH copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.9429734
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.107182
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.219061
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.508868
-----------
Joe copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.772105
-----------
Gumbel copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1
-----------
Gumbel copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.061627
-----------
Gumbel copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.125
-----------
Gumbel copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.285714
-----------
Gumbel copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.428571
-----------
Galambos copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0
-----------
Galambos copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.269205
-----------
Galambos copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.3581411
-----------
Galambos copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.5458864
-----------
Galambos copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.699274
-----------
Husler-Reiss copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
-----------
Husler-Reiss copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.5627681
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
-----------
Husler-Reiss copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.6830113
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
-----------
Husler-Reiss copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.9228086
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
-----------
Husler-Reiss copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 1.111064
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
Warning in log(l.u[, 1L]/u2p) : NaNs produced
Warning in log(u1p/u2p) : NaNs produced
Warning in log(u2p/u1p) : NaNs produced
-----------
Tawn copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0
-----------
Tawn copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.1682746
-----------
Tawn copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.312409
-----------
Tawn copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.5876047
-----------
Tawn copula, dim. d = 2 
Dimension:  2 
Parameters:
  alpha   = 0.7613025
-----------
t-ev copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0
  df      = 4
-----------
t-ev copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.000739897
  df      = 4.000000000
-----------
t-ev copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.2597801
  df      = 4.0000000
-----------
t-ev copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.5545006
  df      = 4.0000000
-----------
t-ev copula, dim. d = 2 
Dimension:  2 
Parameters:
  rho.1   = 0.6818544
  df      = 4.0000000
-----------
> options(op)
> 
> ### FIXME:  Should see nothing below ( <==> length(ccc) == 0) ! ------------------------------
> 
> ccc <- unlist(copChk) # ==> drop all those that have NULL, i.e. *no* error above :
> names(ccc) <- sub(".message$", "", names(ccc))
> structure(as.list(ccc), class="simple.list")
                   _                                
claytonCopula2     is.na(Fu[u.ina]) are not all TRUE
frankCopula2       0 <= Fu[!u.ina] are not all TRUE 
galambosCopula1    is.na(fu[u.iNA]) are not all TRUE
galambosCopula2    0 <= Fu[!u.ina] are not all TRUE 
galambosCopula3    0 <= Fu[!u.ina] are not all TRUE 
galambosCopula4    0 <= Fu[!u.ina] are not all TRUE 
galambosCopula5    0 <= Fu[!u.ina] are not all TRUE 
huslerReissCopula1 0 <= Fu[!u.ina] are not all TRUE 
huslerReissCopula2 0 <= Fu[!u.ina] are not all TRUE 
huslerReissCopula3 0 <= Fu[!u.ina] are not all TRUE 
huslerReissCopula4 0 <= Fu[!u.ina] are not all TRUE 
huslerReissCopula5 0 <= Fu[!u.ina] are not all TRUE 
tawnCopula1        is.na(Fu[u.ina]) are not all TRUE
tawnCopula2        is.na(Fu[u.ina]) are not all TRUE
tawnCopula3        is.na(Fu[u.ina]) are not all TRUE
tawnCopula4        is.na(Fu[u.ina]) are not all TRUE
tawnCopula5        is.na(Fu[u.ina]) are not all TRUE
> 
> stopifnot(length(ccc) <= 17)
> 
> ## claytonCopula2     is.na(Fu[u.ina]) are not all TRUE  ((indepCopula case; progress already))
> ## frankCopula2       0 <= Fu[!u.ina] are not all TRUE   ((indepCopula case; fix?
> 
> ## galambosCopula1    is.na(fu[u.iNA]) are not all TRUE  ((indepCopula case; progress already))
> ## galambosCopula2    0 <= Fu[!u.ina] are not all TRUE
> ## galambosCopula3    0 <= Fu[!u.ina] are not all TRUE
> ## galambosCopula4    0 <= Fu[!u.ina] are not all TRUE
> ## galambosCopula5    0 <= Fu[!u.ina] are not all TRUE
> 
> ## huslerReissCopula1 0 <= Fu[!u.ina] are not all TRUE
> ## huslerReissCopula2 0 <= Fu[!u.ina] are not all TRUE
> ## huslerReissCopula3 0 <= Fu[!u.ina] are not all TRUE
> ## huslerReissCopula4 0 <= Fu[!u.ina] are not all TRUE
> ## huslerReissCopula5 0 <= Fu[!u.ina] are not all TRUE
> 
> ## tawnCopula1        is.na(Fu[u.ina]) are not all TRUE
> ## tawnCopula2        is.na(Fu[u.ina]) are not all TRUE
> ## tawnCopula3        is.na(Fu[u.ina]) are not all TRUE
> ## tawnCopula4        is.na(Fu[u.ina]) are not all TRUE
> ## tawnCopula5        is.na(Fu[u.ina]) are not all TRUE
> 
> ###-------------- "Unit" Tests for specific cases: ------------------------------
> 
> n.c <- as.integer(if(doExtras) 1001 else 101)
> 
> ##== plackettCopula:  density log(dCopula()) vs "better"  dCopula(*, log=TRUE) ----- _UNFINISHED_ TODO
> 
> ##' difference   log(c(..)) -  c(.., log=TRUE)
> mkDfn <- function(u) {
+     stopifnot(is.numeric(u), length(u) == 2 || (is.matrix(u) && ncol(u) == 2))
+     Vectorize(function(alp) { C <- plackettCopula(alp); log(dCopula(u,C)) - dCopula(u,C,log=TRUE) })
+ }
> Df <- mkDfn(u = c(.01,.99))
> curve(Df, 0.99, 1.01, n=n.c) # alpha ~= 1 -- log=TRUE should be better -- but difference seem random!
> curve(Df, 0,    2,    n=n.c) # to my surprise, alpha ~= 0 makes big difference
> curve(Df, 1e-9, .01, log="x")
> curve(abs(Df(x)), 1e-14, .1, log="xy")# and which one is more accurate ??
> # another u[]
> Df <- mkDfn(u = c(.1,.2))
> curve(Df, 0.99, 1.01, n=n.c) # alpha ~= 1 -- log=TRUE: diff. larger for 1+eps than 1-eps
> curve(Df, 0,    2,    n=n.c) # all very small
> 
> # and another
> curve(mkDfn(u = c(.01,.02 ))(x), 0,  2,   n=n.c) # all very small
> curve(mkDfn(u = c(.98,.999))(x), 0,  2,   n=n.c) # growing with alpha
> curve(mkDfn(u = c(.98,.999))(x), .1, 1e4, n=n.c, log="x") # growing with alpha, still |.| < 6e-12
> curve(mkDfn(u = c(.01,.001))(x), .01,1e8, n=n.c, log="x") # (growing) still small: |.| < 6e-15
> 
> ##' 2-dim u[] in the low-density corner
> mku <- function(n, one. = 0.999, eps = c(4^-(2:n2), b.n^-((n-n2):n))) {
+     ## find basis for eps for very small eps
+     stopifnot(n >= 2)
+     n2 <- n %/% 2
+     R.n <- - .Machine$double.min.exp / n
+     b.n <- trunc(100*(2^R.n + 0.05))/100
+     stopifnot(.5 > eps, eps > 0, 0.9 <= one., one. <= 1)
+     cbind(one., eps, deparse.level = 0L)
+ }
> if(doExtras) withAutoprint({
+     alpha <- 1e6 ## <- large alpha
+     (cop <- plackettCopula(alpha))
+     summary(u <- mku(64))
+     d1 <- dCopula(u, cop, log=TRUE)
+     d2 <- log(dCopula(u, cop))
+     summary(d1-d2)
+     all.equal(d1,d2, tol=0)
+ })
> 
> 
> proc.time()
   user  system elapsed 
   4.73    0.20    4.93