Hi Michael,
You may also want to combine the contribution of a group of input variables on 
a group of output variables, eg when the output is the vector of coefficients 
of a Karhunen Loeve expansion... At the end it produces a lot of 
combinations!It may also be interesting to produce the conditional expectations 
as functions as a post-processing of a chaos expansion. Your opinion on this 
point?
Cheers
Régis
    Le mardi 3 juillet 2018 à 17:49:31 UTC+2, BAUDIN Michael 
<[email protected]> a écrit :  
 
 
Hi Régis,
 
  
 
Indeed, my mail is back, but not the data… yet. Hence, I would be happy to show 
a basic Python implementation of the total group Sobol’ indice, but the most 
efficient way is to get a template script … which is in the encrypted hard 
drive of the failed computer.L This should be solved within a couple of days, I 
hope.
 
  
 
Indeed a generic regular expression engine could to the trick, but it is much 
more complicated than required, I think. After all, here is all the choices we 
have :
 
·        First or total indice
 
·        By variable or by group
 
This generates only 4 indices that can be required, to my knowledge :
 
·        first order indice by variable
 
·        total order indice by variable
 
·        first order indice by group
 
·        total order indice by group
 
and this is it.
 
  
 
In attachment, I put three files in order to show a basic test case for what I 
would like to do
 
·        crue-generate-with-dependency.py : generates a sample with 8 variables 
in input and 1 variable in output
 
·        crue-8var-copula.csv
 
·        chaos-by-group.py : creates a chaos from the sample and explores the 
indices from it.
 
This last script ends with the lines :
 
  
 
# Groups of variables (known)
 
allgroups = [[0],[1],[2],[3,4,5],[6,7]]
 
# Create a functional chaos model
 
algo = ot.FunctionalChaosAlgorithm(x, y)
 
print(algo.getAdaptiveStrategy())
 
print(algo.getProjectionStrategy())
 
algo.run()
 
result = algo.getResult()
 
fccr = ot.FunctionalChaosRandomVector(result)
 
for group in allgroups:
 
    s = fccr.getSobolGroupedIndex(group)
 
    # TODO : get the total order indice of each group
 
    print("Group=%s, S=%.2f" % (str(group),s))
 
  
 
  
 
The script I have in my hard drive allow to loop over the multi-indices, 
computing the required indices within the loop : I will complete my message 
when possible. Meanwhile, I think that there is an issue within the API since 
it does not seem possible to get the enumeration rule from the functional chaos 
“algo” : am I wrong ?
 
  
 
Your analysis of the estimate of these grouped indices by sampling make it 
clear that this is an advantage of the chaos : estimating the group indices is 
much easier with chaos than by sampling. 
 
  
 
Thank you very much with the copula decomposition script : I think that I can 
use almost as is !
 
  
 
Best regards,
 
  
 
Michaël
 
  
 
  
 
De : [email protected] 
[mailto:[email protected]]
Envoyé : mardi 19 juin 2018 10:40
À : [email protected]; BAUDIN Michael <[email protected]>
Objet : Re: [ot-users] Grouped Sobol indices
 
  
 
Hi Michael,
 
  
 
Nice to see that your email is back ;-)
 
  
 
Concerning your first point and a part of your post-scriptum, an idea could be 
to provide a kind of regular expression based filter to compute the indices you 
want. Given a set of input variables, you may want to collect all the 
coefficients involving each of the variables within this set with a positive 
exponent, or all the coefficients involving at least one of these variables 
with a positive exponent, or at least all these variables with a positive 
exponent, or all the coefficients involving none of these variables with a 
positive exponent and so on. And you may want to collect these coefficients 
output component by output component, or for a given set of output component, 
or... You see what I mean? At one point we have to bound the number of methods, 
so a generic extraction mechanism and a dedicated set of filters may be the 
most efficient way to go. Do you have a list of indices of interest, with their 
formal definition in terms of conditional expectation? It would help to see how 
to implement this method.
 
  
 
Concerning the estimation of group indices by sampling, as soon as you are able 
to express them in terms of variance of conditional expectations, it is quite 
straightforward to implement a scheme similar to the Saltelli method. The main 
difficulties are the potentially large number of indices to compute, and the 
huge number of simulations required to sample conditional expectations with a 
large conditioning dimension.
 
  
 
For the last point, there is currently no way to inspect a given multivariate 
distribution to see how to factor it into smaller pieces. The point of view in 
OpenTURNS is the opposite: you can build such a product multivariate 
distribution using the ComposedCopula class (each sub-copula describes the 
dependence within a group of dependent variables), but there is no easy way to 
perform the inverse computation. From a C++ perspective, you can ask for the 
copula of your distribution and try to dynamic_cast it into a ComposedCopula to 
get the associated ComposedCopula, but if it is the case, odds are large that 
you could have access to this information by looking at your Python script a 
few lines above!
 
Nevertheless, if you really need this kind of service, it can be implemented 
quite easily using Csiszar divergences. In the following script, I use the 
(squared) Hellinger divergence to check if a given distribution is a product of 
two lower dimensional distributions at a given position. The computation of the 
divergence is done by sampling (Monte Carlo). If you have *exact* independence, 
the sampling size can be very low (less than 100) as the divergence is 
*exactly* zero and the kernel function drops to nearly zero at each point of 
the sample. From this script, it is straightforward to get the associated block 
distributions.
 
  
 
import openturns as ot
 
import numpy as np
 
from time import time
 
  
 
def isProduct(distribution, cut, N = 10000, epsilon = 1.0e-16):
 
    dim = distribution.getDimension()
 
    distribution1 = distribution.getMarginal([i for i in range(cut)])
 
    distribution2 = distribution.getMarginal([cut+i for i in range(dim-cut)])
 
    def kernel(x):
 
        x = np.array(x)
 
        s = x.shape
 
        x1 = x[:, 0:cut]
 
        pdf1 = np.array(distribution1.computePDF(x1))
 
        x2 = x[:, cut:dim]
 
        pdf2 = np.array(distribution2.computePDF(x2))
 
        pdf = np.array(distribution.computePDF(x))
 
        return ((1.0 - np.sqrt(pdf1 * pdf2 / pdf))**2 * pdf).reshape(s[0], 1)
 
  
 
    value = ot.PythonFunction(dim, 1, 
func_sample=kernel)(distribution.getSample(N)).computeMean()[0]
 
    return (value < epsilon)
 
  
 
def decompose(distribution, N = 10000, epsilon = 1.0e-16):
 
    dim = distribution.getDimension()
 
    blocks = list()
 
    currentBlock = [0]
 
    for i in range(1, dim-1):
 
        if isProduct(distribution, i, N, epsilon):
 
            blocks.append(currentBlock)
 
            currentBlock = [i]
 
        else:
 
            currentBlock.append(i)
 
    currentBlock.append(dim-1)
 
    blocks.append(currentBlock)
 
    return blocks
 
  
 
if __name__ == "__main__":
 
    d1 = 2
 
    d2 = 1
 
    d3 = 4
 
    dimension = d1 + d2 + d3
 
    R = ot.CorrelationMatrix(dimension)
 
    for i in range(d1):
 
        for j in range(i):
 
            R[i, j] = 0.1 / (1 + i + j)
 
    for i in range(d2):
 
        for j in range(i):
 
            R[i + d1, j + d1] = 0.1 / (1 + i + j)
 
    for i in range(d3):
 
        for j in range(i):
 
            R[i + d1 + d2, j + d1 + d2] = 0.1 / (1 + i + j)
 
    print(R)
 
  
 
    t0 = time()
 
    print("Test 1: Normal distribution")
 
    distribution = ot.Normal([0.0]*dimension, [1.0]*dimension, R)
 
    print("blocks=", decompose(distribution))
 
    print("Test 2: Composed copula")
 
    distribution = ot.ComposedCopula([ot.ClaytonCopula(5.0)]*5 + 
[ot.NormalCopula(R)])
 
    print("blocks=", decompose(distribution))          
 
    print("t=", time() - t0, "s")
 
  
 
  
 
Cheers
 
  
 
Régis LEBRUN
 
  
 
  
 
  
 
  
 
Le lundi 18 juin 2018 à 16:53:20 UTC+2, BAUDIN Michael <[email protected]> 
a écrit : 
 
  
 
  
 
Hi !
 
 
 
I have a computer code with dependent input vector and I would like to perform 
sensitivity analysis on the vector output. To do this, I would like to estimate 
the Sobol’ indices. In order to workaround the dependent inputs, I would like 
to gather dependent input variables into independent groups and to estimate the 
Sobol’ indices on these groups.
 
 
 
·        My first idea is to use a chaos expansion in order to estimate Sobol’ 
indices. It is easy to estimate the first order Sobol’ indices with the 
getSobolGroupedIndex method of the FunctionnalChaosSobolIndices class:
 
 
 
http://openturns.github.io/openturns/master/user_manual/response_surface/_generated/openturns.FunctionalChaosSobolIndices.html
 
 
 
However, the corresponding total order indices method seem to be unavailable: 
there is no getSobolGroupedTotalIndex ?
 
 
 
If so, this is a pity, since this is just another “simple” arrangement of the 
coefficients from the chaos expansion. It would allow to estimate the total 
effect of a group. Compared to the group Sobol first indice, it would allow to 
see if there are interactions with this group and other variables outside from 
the group. Also, if the group Sobol total indice is close to zero, this means 
that this group has no effect of the variability of the output. 
 
 
 
·        The second topic of my post is the SobolIndicesAlgorithm class :
 
 
 
http://openturns.github.io/openturns/master/user_manual/_generated/openturns.SobolIndicesAlgorithm.html
 
 
 
I can see that there is no “group” method here. Is there a theoretical way of 
computing grouped Sobol’ indices with these estimators? 
 
 
 
·        Finally, I would like to know if there is a way of using the joint 
distribution in order to generate the independent groups of dependent input 
variables? In other words, with a given multivariate joint distribution, is 
there a way of generating a list of sub-lists of indices such that:
 
o  Each sub-list contain dependent variables,
 
o  Two sub-lists are independent.
 
With such a goal, what classes / methods to use ?
 
 
 
Best regards,
 
 
 
Michaël
 
 
 
PS
 
By the way, the total Sobol’ indices of a group of variables lead to a way of 
removing input variables from the model. Indeed, we could compute the largest 
group having a total Sobol’ indice lower than a given small threshold (says 
0.01 for example). By design, all variables from this group could be replaced 
by fixed inputs without changing the variability of the output. Furthermore, if 
we had the distribution of this estimator, we could guarantee this property 
with, say, 95% confidence (i.e. including the variability from the estimate).
 
In order to compute this group, the variables could first be ordered by 
decreasing total Sobol’ indices. Then, we could add the variables into the 
group until the total Sobol’ indice exceeds the threshold: the final group is 
the one just before this step.
 
With this raw algorithm, the chaos expansion has a great advantage: given that 
the decomposition is known, the estimate is almost CPU-free. 
 

Ce message et toutes les pièces jointes (ci-après le 'Message') sont établis à 
l'intention exclusive des destinataires et les informations qui y figurent sont 
strictement confidentielles. Toute utilisation de ce Message non conforme à sa 
destination, toute diffusion ou toute publication totale ou partielle, est 
interdite sauf autorisation expresse.
 
Si vous n'êtes pas le destinataire de ce Message, il vous est interdit de le 
copier, de le faire suivre, de le divulguer ou d'en utiliser tout ou partie. Si 
vous avez reçu ce Message par erreur, merci de le supprimer de votre système, 
ainsi que toutes ses copies, et de n'en garder aucune trace sur quelque support 
que ce soit. Nous vous remercions également d'en avertir immédiatement 
l'expéditeur par retour du message.
 
Il est impossible de garantir que les communications par messagerie 
électronique arrivent en temps utile, sont sécurisées ou dénuées de toute 
erreur ou virus.
____________________________________________________
 
This message and any attachments (the 'Message') are intended solely for the 
addressees. The information contained in this Message is confidential. Any use 
of information contained in this Message not in accord with its purpose, any 
dissemination or disclosure, either whole or partial, is prohibited except 
formal approval.
 
If you are not the addressee, you may not copy, forward, disclose or use any 
part of it. If you have received this message in error, please delete it and 
all copies from your system and notify the sender immediately by return message.
 
E-mail communication cannot be guaranteed to be timely secure, error or 
virus-free.
 
_______________________________________________
OpenTURNS users mailing list
[email protected]
http://openturns.org/mailman/listinfo/users
 



Ce message et toutes les pièces jointes (ci-après le 'Message') sont établis à 
l'intention exclusive des destinataires et les informations qui y figurent sont 
strictement confidentielles. Toute utilisation de ce Message non conforme à sa 
destination, toute diffusion ou toute publication totale ou partielle, est 
interdite sauf autorisation expresse.

Si vous n'êtes pas le destinataire de ce Message, il vous est interdit de le 
copier, de le faire suivre, de le divulguer ou d'en utiliser tout ou partie. Si 
vous avez reçu ce Message par erreur, merci de le supprimer de votre système, 
ainsi que toutes ses copies, et de n'en garder aucune trace sur quelque support 
que ce soit. Nous vous remercions également d'en avertir immédiatement 
l'expéditeur par retour du message.

Il est impossible de garantir que les communications par messagerie 
électronique arrivent en temps utile, sont sécurisées ou dénuées de toute 
erreur ou virus.
____________________________________________________

This message and any attachments (the 'Message') are intended solely for the 
addressees. The information contained in this Message is confidential. Any use 
of information contained in this Message not in accord with its purpose, any 
dissemination or disclosure, either whole or partial, is prohibited except 
formal approval.

If you are not the addressee, you may not copy, forward, disclose or use any 
part of it. If you have received this message in error, please delete it and 
all copies from your system and notify the sender immediately by return message.

E-mail communication cannot be guaranteed to be timely secure, error or 
virus-free.
  
_______________________________________________
OpenTURNS users mailing list
[email protected]
http://openturns.org/mailman/listinfo/users

Reply via email to