Hi Michael,
Have a look 
here:https://github.com/openturns/www/blob/master/_images/scripts/plot_kriging.py

I will check your script next rainy day!
Cheers
Régis


    Le vendredi 10 août 2018 à 14:22:58 UTC+2, BAUDIN Michael 
<[email protected]> a écrit :  
 
  
Hi !
 
  
 
I am currently working on kriging and try to reproduce the graphics which is 
printed in the carrousel of openturns.org :
 
  
 
https://www.qwant.com/?q=openturns%20kriging&t=images&o=0:1e0b94f86add70a402dfc53f500aa057&size=all&license=all&freshness=all&color=all&imagetype=all&source=web
 
  
 
More precisely, I would like to compute the 95% confidence bounds of the 
kriging predictions. However, I could not find the Python script corresponding 
to the figure neither on openturns.org nor in the doc (nor in the Slides of the 
User’s days, nor in the Use Cases Guide). Is this Python script available ?
 
  
 
I also tried to produce the script by myself, and failed.... In PS, I put a 
Python script which tries to reproduce the figure. The confidence bounds 
produced by the script are extremely tight compared to what I expected. My 
guess is that the issue is related to the estimate of the covariance model : in 
OpenTURNS, the SquaredExponential is anisotropic while the figure might be 
produced with an isotropic covariance. I also tried to prevent the optimization 
algorithm from tuning the theta parameter with the setOptimizeParameters 
method, but I wasn’t able to use correctly (the parameters still change, 
however with a different value). What is the explanation for these tight bounds 
?
 
  
 
Best regards,
 
  
 
Michaël
 
PS
 
  
 
# coding: utf-8
 
#
 
# # Références
 
#
 
# Metamodeling with Gaussian processes, Bertrand Iooss, EDF R&D, 2014, 
www.gdr-mascotnum.fr/media/sssamo14_iooss.pdf
 
  
 
#
 
import numpy as np
 
import openturns as ot
 
import pylab as pl
 
  
 
#
 
n_test = 100
 
x_test_coord = np.linspace(0,12,n_test)
 
x_test = ot.Sample(x_test_coord,1)
 
y_test = np.sin(x_test)
 
  
 
# Base d'apprentissage de 7 points.
 
#
 
x_train_coord = [1.,3.,4.,6.,7.9,11., 11.5]
 
x_train = ot.Sample(x_train_coord,1)
 
y_train = np.sin(x_train)
 
n_train = len(x_train)
 
n_train
 
  
 
#
 
dimension = 1
 
basis = ot.ConstantBasisFactory(dimension).build()
 
covarianceModel = ot.SquaredExponential([1.]*dimension, [5.0])
 
algo = ot.KrigingAlgorithm(x_train, y_train, covarianceModel, basis)
 
algo.setOptimizeParameters(False)
 
algo.run()
 
result = algo.getResult()
 
krigeageMM = result.getMetaModel()
 
  
 
#
 
krigeageMM = result.getMetaModel()
 
y_test_MM = krigeageMM(x_test)
 
  
 
#
 
level = 0.95
 
alpha = (1-level)/2
 
sampleQinf = ot.Sample(n_test,1)
 
sampleQsup = ot.Sample(n_test,1)
 
for i in range(n_test):
 
    xi = x_test[i]
 
    var = result.getConditionalCovariance(xi)[0,0]
 
    yi = y_test_MM[i][0]
 
   N = ot.Normal(yi,np.sqrt(var))
 
    q025 = N.computeQuantile(alpha)[0]
 
    sampleQinf[i,0] = q025
 
    q975 = N.computeQuantile(1-alpha)[0]
 
   sampleQsup[i,0] = q975
 
  
 
#
 
pl.plot(x_test,y_test,"k-",label="Exact")
 
pl.plot(x_test,y_test_MM,"b--",label="Krigeage")
 
pl.plot(x_train,y_train,"ro",label="Apprentissage")
 
pl.plot(x_test,sampleQinf,"g-",label="%.2f%%" % (100*level))
 
pl.plot(x_test,sampleQsup,"g-")
 
pl.title("sin")
 
pl.xlabel("X")
 
pl.ylabel("Y")
 
pl.legend()
 
  
 
  
 



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
  
_______________________________________________
OpenTURNS users mailing list
[email protected]
http://openturns.org/mailman/listinfo/users

Reply via email to