Vinecop pdf

[tnagler]

While implementing the pdf I noticed that values in low-density regions were nan. This was caused by rounding errors when evaluating the pair-copula pdf and h-hfunctions too close to the boundaries. So it seems as if the truncation in VineCopula is necessary.

I modified rotate_u and renamed it cut_and_rotate. It cuts off all values outside the interval [1e-20, 1 -1e20] (VineCopula uses 1e-10). This function is now also called when rotation = 0.

[tvatter]

1. In my opinion, cut_and_rotate should be two separate functions. Furthermore, can you benchmark the overhead caused by using cut_and_rotate when rotation_ = 0 in the estimation routine? Calling cut_and_rotate when rotation_ = 0 means that we're copying the data every-time we compute the likelihood, which is not very wise when optimizing it.
2. As explained on Skype, I think that the pdf implementation for vinecop objects does not look very "C++ish" (i.e., I understand that it is a direct translation of the C code). I think that we should really try to find an alternative way of computing this, using "vectorized" operations. I would still merge the PR for now (after answers to my first comment above), but we definitely need to think about it.

[tnagler]

1. I thought of that yesterday too. Two separate functions is certainly neater. But I don't think there's a way around copying when we want to cut. We could modify the original data, but this is not really a clean solution.

2. We could vectorize over i, but this will exceed the memory in very high dimensions. I "unvectorized" the VineCopula version exactly for that reason in the December release.

[tvatter]

1. I would still be interested in a benchmark tough. An alternative would be to simply avoid calling the Bicop::loglik in ParBicop::fit. We could either use computing the log-likelihood "manually" in mle_objective and pmle_objective or do something similar as for the pdf:

double Bicop::loglik(const MatXd& u)
{
    return loglik_default(cut_and_rotate(u));
}
double Bicop::loglik_default(const MatXd& u)
{
    VecXd ll = this->pdf_default(u);
    ll = ll.array().log();
    return ll.sum();
}

I have a preference for the later (to be consistent with pdf), what do you think?
2. I implemented the vectorized version in commit 27f235f.

[tnagler]

1. I would not call Bicop::loglik in ParBicop::fit but rather cut once and then call this->pdf_default(u).array().log().sum(). If we don't need it anywhere else there's no point in having an own function for loglik_default.

2. Perfect!

[tvatter]

1. This won't work since pdf_default is a private method. In fact, my two solutions won't work either for the same reason... Anyway, I did benchmark Bicop::fit when copying and without copying, and the overhead seems very low. Let's keep it that way for now and think about this as a potential area of improvement later!