Utiliser GridGain pour le calcul de la Value At Risk

Après un premier article introduisant l’intérêt de la Value At Risk and du calcul en grille, nous allons désormais étudier l’implémentation de cet algorithme en utilisant un middleware de calcul en grille. J’ai choisi GridGain, un middleware open source qui implémente le pattern map/reduce (cf. mon précédent article). Pour commencer, je vais donner un aperçu de l’implémentation de la Value At Risk indépendamment de l’architecture de calcul en grille. Ensuite, je décrirai le middleware GridGain et les classes à implémenter pour tirer parti de la grille. Enfin, je fournirai quelques mesures de performances prises sur mon portable et proposerai quelques pistes d’optimisation.

Certaines portions du code ont été supprimées pour des raisons de lisibilité et remplacées par des commentaires (//...).
Le diagramme de classes ci-dessous donne un aperçu des différentes classes :

OptionPricer, Parameters, ParametersValueGenerator et VarComputerMaanger sont destinées au calcul de la Value At Risk. VarComputerGridRunner et VarComputerGridTask sont les classes implémentées pour tirer parti de la grille. Configuration est partagée à travers toutes les couches.

Implémentation de la Value at Risk

Cette classe calcule le prix de l’option en fonction des paramètres fournis. Il s’agit strictement de l’implémentation de la formule de Black and Scholes. Les fonctions mathématiques comme le calcul d’une distribution normale (NormalDistribution) sont fournies par la librairie Apache commons Math.

public class OptionPricer implements Serializable {
	private static final long serialVersionUID = 5294557359966059686L;
	private static int NUMBER_OF_TRADE_DAY = 252;
	
	public static class Parameters implements Serializable {
		/**  For Serialization	 */
		private static final long serialVersionUID = -3967053376944712365L;
		/** Time in days before the maturity	 */
		private int t;
		/** Time in year before the maturity  */
		private double tInYear;
		/** Spot (actual price) of the underlying  */
		private double s0;
		/**  Strike (target price) of the underlying  */
		private double k;
		/**  Risk free interest rate  */
		private double r;
		/**  Implicit volatility of the Equity	 */
		private double sigma;
		
		public Parameters(int t, double s0, double k, double r, double sigma) {//...	}
		public Parameters(final Parameters parameters) { //...	}
		
		//Getters and Setters...
	}//End Parameters class
	
	 /** Initialize only one time and reuse the distribution	 */
	private NormalDistribution normalDistribution = new NormalDistributionImpl(0, 1);
	private final Parameters initialParameters;
	private Parameters parameters;
	
	public OptionPricer(final Parameters intialParameters) {
		this.initialParameters = intialParameters;
		this.parameters = new Parameters(this.initialParameters);
	}
	//Getter and Setter both for initialParameters and parameters...

	/**
	 * d_1 = \frac{1}{\sigma\sqrt{t}} \left[ \ln \left( \frac{S_0}{K} \right) + \left( r + \frac{1}{2}\sigma^2 \right)t \right]
	 * 
	 * @param spot actualized by a the yield for this simulation
	 */
	private double computeD1(final double s0) {
		//denominator
		double denom = parameters.sigma * Math.sqrt(parameters.getTInYear());
		//numerator
		//I don't explain - 0.5 * sigma^2 but in the Excel sheet it is like that
		double num = Math.log(s0 / parameters.k) + (parameters.r - 0.5 * parameters.sigma * parameters.sigma) * parameters.getTInYear(); 
		return num / denom;
	}
	
	/**
	 * d_2 = d_1 - \sigma \sqrt{t}
	 * @param d1
	 * @return
	 */
	private double computeD2(final double d1) {
		double result = d1 - parameters.sigma * Math.sqrt(parameters.getTInYear()); 
		return result;
	}

	/**
	 * Compute the price of the call
	 * C(S_0,K,r,t,\sigma) = S_0 \mathcal{N}(d_1) - K e^{-rt}\mathcal{N}(d_2)
	 * \mathcal{N}(x) = \int_{-\infty}^{x} \frac{1}{\sqrt{2\pi}}e^{-\frac{1}{2}u^2} du
	 * 
	 * @return
	 * @throws MathException 
	 * 
	 * @param spot yield Allows to build a actual yield for each simulation
	 */
	public double computeCall() throws MathException {
		final double s0 = parameters.s0;
		final double d1 = computeD1(s0);
		final double d2 = computeD2(d1);
		double result = s0 * normalDistribution.cumulativeProbability(d1) - parameters.k * Math.exp(-1 * parameters.r * parameters.getTInYear()) * normalDistribution.cumulativeProbability(d2);
		return result;
	}
}

Deux points méritent d’être signalés. Tout d’abord, pour des raisons de performance, j’ai choisi d’utiliser des nombres à virgules flottantes. De façon à éviter certains écueils, certaines précautions doivent être prises. Aucun arrondi ou fonction non linéaire ne doit avoir lieu durant le calcul. Toutes les divisions doivent être reportées à la fin du calcul afin de ne pas diminuer la précision du calcul. Le second point notable est la définition, toujours pour des raisons de performance, d’un champ initialParameters et d’une référence au champ parameter. Je détaillerai ce point ci-dessous en m’aidant de la classe suivante.

public class ParametersValueGenerator implements Serializable {
	private static final long serialVersionUID = 2923455335542820656L;
	private final double underlyingYieldHistoricalVolatility;
	private final double interestRateVariationVolatility;
	private final double implicitVolatilityVolatility;
	
	/** Used to generate the raws */
	private RandomData randomGenerator = new RandomDataImpl();
	/** Reference to OptionPricer initialParameters */
	private OptionPricer.Parameters initialParameters;

	/**
	 * @param underlyingYieldHistoricalVolatility
	 *            Volatity of the yield of the undelying price (it is the
	 *            standard deviation for the normal distribution)
	 * @param interestRateVariationVolatility
	 *            Volatility of the risk free interest rate (it is the
	 *            standard deviation for the normal distribution)
	 * @param implicitVolatilityVolatility
	 *            Volatility of the underlying volatility (we probably take
	 *            the standard deviation of an index like VIX
	 * @see http://en.wikipedia.org/wiki/VIX
	 */
	public ParametersValueGenerator(
			double underlyingYieldHistoricalVolatility,
			double interestRateVariationVolatility,
			double implicitVolatilityVolatility) {
		this.underlyingYieldHistoricalVolatility = underlyingYieldHistoricalVolatility;
		this.interestRateVariationVolatility = interestRateVariationVolatility;
		this.implicitVolatilityVolatility = implicitVolatilityVolatility;
	}

        /**
	    Reseed with a value that is more likely to be unique
		This is mandatory to avoid generating two identical random number series
		To be really sure to not have the same value in two different thread or computer, 
		I use a secure random generator instead of the current millisecond
		
		This method is time consumming, use it a minimum number of times
	 */
	public void intializeRandomGenerator() {
		randomGenerator.reSeed(securedRandomGenerator.nextLong());
	}
	
	/**
	 * Set the reference
	 * @param initialParameters
	 * @param parameters
	 */
	public void initializeReference(final OptionPricer.Parameters initialParameters) {
		this.initialParameters = initialParameters;
	}

	/**
	 * Put random generated values into parameters instance by modifying it
	 * by reference
	 * 
	 * @param parameters
	 */
	public void setNewGeneratedParameters(OptionPricer.Parameters parameters) {
		if(parameters == null) { throw new IllegalArgumentException("parameters should not be null"); }
		if(this.initialParameters == null) { throw new IllegalStateException("this.initialParameters should not be null"); }
		
		if(this.underlyingYieldHistoricalVolatility != 0) {
			parameters.setS0(this.initialParameters.getS0() * (1+randomGenerator
				.nextGaussian(0, this.underlyingYieldHistoricalVolatility)));
		}
		if(this.interestRateVariationVolatility != 0) {
			parameters.setR(this.initialParameters.getR() + randomGenerator
					.nextGaussian(0, this.interestRateVariationVolatility));
		}
		if(this.implicitVolatilityVolatility != 0) {
			parameters.setSigma(randomGenerator
					.nextGaussian(this.initialParameters.getSigma(),
							this.implicitVolatilityVolatility));
		}
	}
	//toString() method...
}

Cette classe génère des valeurs aléatoires pour la méthode de Monte-Carlo. Cette implémentation élémentaire génère des données suivant une distribution normale. Cependant, on notera d’un point de vue financier, que la distribution normale est utilisée pour modéliser non pas directement les valeurs mais certaines de leurs caractéristiques :

le rendement du prix du sous-jacent (conduisant à value=initial*(1+nextGaussian()))
la variation du taux d’intérêt (value=initial+nextGaussian())
la valeur de la volatilité implicite (value=nextGaussian())

La méthode setNewGeneratedParameters() modifie par référence le champ parameters de l’instance OptionPricer en utilisant pour son calcul une référence à initialParameters. La méthode initializeReference() permet de la faire pointer sur l’objet initialParameters de l’OptionPricer. Il aurait été plus simple de renvoyer un nouvel objet Parameters pour chaque simulation, mais quelques mesures avec JVisualVM montrent que le constructeur de la classe Parameters était l’une des méthodes les plus consommatrices. En effet, allouer la mémoire sur le tas nécessite de rechercher le prochain emplacement avec suffisamment de place. C’est une opération coûteuse. Passer chaque valeur native à la méthode de l’OptionPricer aurait potentiellement été plus efficace. Mais ce passage par référence semblait être un bon compromis entre la performance et la maintenabilité.

Un autre point à mentionner est la méthode intializeRandomGenerator () . Sans elle, vous obtiendrez un comportement difficilement explicable qui peut vous laisser perplexe pendant plusieurs heures. Je vous décrirai le problème et la solution un peu plus tard, mais pour cela je dois d’abord rentrer un peu plus dans les détails d’implémentation.

L’outil GridGain et l’implémentation Map/Reduce

L’outil fournit un certain nombre de fonctionnalités permettant d’exécuter le code sur une grille sans avoir à coder la logique bas niveau de distribution. Parmi les fonctionnalités fournies « out of the box », nous pouvons citer :

La découverte automatique des noeuds disponibles
Le déploiement du code et des données sur chaque noeud
La répartition de charge (ou load balancing)
La reprise sur erreur (ou failover) si un noeud est en erreur ou ne répond plus

Concrètement, la configuration par défaut de GridGain se démarre très simplement dans un environnement de développement. Le seul pré-requis est que GridGain soit déployé et exécuté sur chaque noeud. L’application conçue démarre un noeud embarqué. Après la découverte des autres noeuds, la grille est initialisée comme le montre le contenu de la console ci-dessous (dans mon cas une grille avec 4 noeuds sur une machine équipée d’un processeur quad-core):

>>> -------------------
>>> Discovery Snapshot.
>>> -------------------
>>> Number of nodes: 4
>>> Topology hash: 0xAEFB9485
>>> Local: CD766E36-2A09-45B5-895F-860F5DFBA386, 10.1.106.135, Windows 7 x86 6.1
, mbo, Java(TM) SE Runtime Environment 1.6.0_20-b02
>>> Remote: B5278A02-CD63-4E7E-89FF-BCA791E9045E, 10.1.106.135, Windows 7 x86 6.
1, mbo, Java(TM) SE Runtime Environment 1.6.0_20-b02
>>> Remote: A4201237-6838-44A3-8B32-3E8B5B1F2FFF, 10.1.106.135, Windows 7 x86 6.
1, mbo, Java(TM) SE Runtime Environment 1.6.0_20-b02
>>> Remote: C9BCB5C8-3C87-4106-9E89-8ABF11FEC959, 10.1.106.135, Windows 7 x86 6.
1, mbo, Java(TM) SE Runtime Environment 1.6.0_20-b02
>>> Total number of CPUs: 4

Le code est déployé automatiquement sur chaque noeud et exécuté.

Le pattern map/reduce

GridGain implémente le pattern Map/Reduce qui est nécessaire pour pouvoir distribuer le calcul. En quelques mots, ce dernier nous vient de la programmation fonctionnelle. Le calcul est découpé en une succession d’un même petit calcul élémentaire appliqué à un grand nombre de données (fonction map). Chaque petit calcul peut ainsi être distribué. Ensuite, les résultats sont fusionnés (fonction reduce). Une description plus approfondie avec un exemple est fournie dans mon précédent article.

L’implémentation du pattern map/reduce

Dans ce cas particulier, le code est organisé de cette façon.
– La classe VarComputerManager implémente la logique de calcul de la VAR dans la méthode computeVar(). Elle n’est pas liée à la grille et implémente la logique métier. Pour être exact, cette méthode ne renvoie par une seule valeur (LA VAR) mais tous les prix inférieurs à la VAR. Je vous renvoie vers mon précédent article pour plus de précision sur la signification métier et l’algorithme de calcul de la VAR.

public class VarComputerManager {
	/**
	 * Compute the VAR of the option defined in the option pricer
	 * @param optionPricer Instance should be instanciated with the correct parameters
	 * @param valueGenerator Instance should be instanciated with the volatilities, and the reference to the default parameters should have been set 
	 * @param drawsNb Number of draws produced on that node
	 * @param totalDrawsNb Total number of draws made on the whole grid
	 * @param varPrecision Between 0 and 1. 0.99 gives the VAR (the percentile) to 1%
	 * @return The n lowest prices generated on that node where n = (1-varPrecision)*totalDrawsNb
	 * @throws MathException
	 */
	public SortedSet<double> computeVar(final OptionPricer optionPricer,
			final ParametersValueGenerator valueGenerator, final int drawsNb, final int totalDrawsNb,
			final double varPrecision, final Configuration configuration) throws MathException {
		final double nbValInPercentile = (totalDrawsNb * (1 - varPrecision));
		final SortedSet</double><double> smallerPrices = new TreeSet</double><double>();
		
		valueGenerator.intializeRandomGenerator();

		for (int i = 0; i < drawsNb; i++) {
			valueGenerator.setNewGeneratedParameters(optionPricer.getParameters());
			final double price = optionPricer.computeCall();

			// For each draw, put the price in the sorted set
			smallerPrices.add(price);
			if(configuration.isCombineAllowed()) {
				// If the set exceeds nbValInPercentile, remove the highest value
				if (smallerPrices.size() > nbValInPercentile) {
					smallerPrices.remove(smallerPrices.last());
				}
			}
		}

		return smallerPrices;
	}
	
	/**
	 * 
	 * @param smallerPrices Modified by reference
	 * @param drawsNb
	 * @param varPrecision
	 */
	public void extractPercentile(SortedSet</double><double> smallerPrices, final int drawsNb, final double varPrecision) {
		final double nbValInPercentile = (drawsNb * (1 - varPrecision));
		// If the set exceeds nbValInPercentile, remove the highest value
		while (smallerPrices.size() > nbValInPercentile) {
			smallerPrices.remove(smallerPrices.last());
		}
	}
}
</double>

Un nombre fixé de tirages de variable aléatoire est exécuté par le code valueGenerator.setNewGeneratedParameters(optionPricer.getParameters()); et le prix du call correspondant aux paramètres de ce tirage est calculé. Les valeurs les plus faibles sont identifiées.
Cet exemple utilise simplement un SortedSet pour les identifier. Un effet de bord de cette méthode est que, dans le cas de deux valeurs de VAR identiques générées par des combinaisons de paramètres différentes, la seconde valeur sera ignorée (add() est une opération optionnelle). En pratique, j’ai constaté ce comportement dans 0,02% des tirages, conduisant à une valeur de VAR potentiellement inférieure à la réalité. Dans le cadre de cet exemple, j’ai choisi de conserver cette imprécision.

– Ensuite, la classe VarComputerGridTask implémente la logique technique pour déployer le code sur la grille. La fonction split() est un synonyme pour le map() du pattern map/reduce. Elle définit comme découper la fonction computeVar() en plus petites fonctions.

public class VarComputerGridTask extends GridTaskSplitAdapter<gridifyargument , SortedSet<Double>> {
	private static final long serialVersionUID = 2999339579067614574L;
	private static String SESSION_CONF_KEY = "SESSION_CONF_KEY";
	private static Logger log = LoggerFactory.getLogger(VarComputerGridTask.class);
	@GridTaskSessionResource
    private GridTaskSession session = null;

	@Override
	protected Collection< ? extends GridJob> split(int gridSize, final GridifyArgument arg)
			throws GridException {
		// Split number of iterations.
        Integer iterations = ((Integer)arg.getMethodParameters()[2]);
       
        // Note that for the purpose of this example we perform a simple homogeneous
        // (non weighted) split on each CPU/Computer assuming that all computing resources 
        // in this split will be identical. 
        int slotsNb = 0;
        //Retrieve the session parameter
        Object[] params = arg.getMethodParameters();
        Configuration conf = (Configuration)params[5];
        //Put it into session to retrieve it in reduce task
        session.setAttribute(SESSION_CONF_KEY, conf);
        if(conf.isThereSeveralSlotsByProcess()) {
        	Grid grid = GridFactory.getGrid();
            for(GridNode node : grid.getAllNodes()) {
            	log.debug("Node {} detected with {} cores", node.getPhysicalAddress(), node.getMetrics().getAvailableProcessors());
            	//1 process per node is presumed. If several processes are running on the machine, the number of cores will be over estimated
            	slotsNb += node.getMetrics().getAvailableProcessors();
            }
            log.debug("{} slots available", slotsNb);
        }
        else {
        	slotsNb = gridSize;
        	log.debug("1 slot by process so {} slots", gridSize);
        }

        // Number of iterations should be done by each slot.
        int iterPerNode = Math.round(iterations / (float)slotsNb);

        // Number of iterations for the last/the only slot.
        int lastNodeIter = iterations - (slotsNb - 1) * iterPerNode;

        List<gridjobadapter <Integer>> jobs = new ArrayList</gridjobadapter><gridjobadapter <Integer>>(gridSize);
        
        for (int i = 0; i < slotsNb; i++) {
            // Add new job reference to the split.
            jobs.add(new GridJobAdapter<Integer>(i == slotsNb - 1 ? lastNodeIter : iterPerNode) {
				private static final long serialVersionUID = 1L;

				/*
                 * Executes grid-enabled method with passed in argument.
                 */
                public TreeSet<double> execute() throws GridException {
                    Object[] params = arg.getMethodParameters();
                    TreeSet</double><double> result = null;                    

                    VarComputerManager varComputerManager = new VarComputerManager();
                    try {
                    	result = (TreeSet</double><double>)varComputerManager.computeVar((OptionPricer)params[0],
								(ParametersValueGenerator)params[1], 
								getArgument(), //Give the split size on that node
								(Integer)params[3],//Give the total size on the grid
								(Double)params[4],
								(Configuration)params[5]);
					} catch (MathException e) {
						log.error("Error computing option pricing", e);
						result = new TreeSet</double><double>();
					}
					
					log.debug("{}", String.format("SlotId %d: Map %d elements on %d and produce %d results", this.hashCode(), getArgument(), (Integer)params[2], result.size()));
					return result;
                }
            });//End GridJobAdapter
        }//End for
        
		return jobs;
	}
</double></gridjobadapter></gridifyargument>

La méthode split() reçoit en argument la taille de la grille et un objet GridArgument contenant les arguments de la méthode VarComputerManager.computeVar() vue ci-dessus. La fonction détermine le nombre de créneaux disponibles selon le nombre de processeurs. Comme nos tâches sont fortement consommatrices de calcul, la répartition la plus efficace est d’un créneau par coeur. Prenons l’exemple de 1000 tirages pouvant être répartis en 10 tâches appliquant chacune la méthode computeVar() sur 100 tirages. Pour GridGain, une tâche est représentée par une instance GridJobAdapter .
De façon à être capable de mesurer l’influence de différentes stratégies (par exemple un thread par processus mais plusieurs processus ou un processus avec plusieurs threads par machine), certaines caractéristiques ont été définies dans un objet Configuration. La GridTaskSession permet de la rendre accessible à chaque méthode map() et reduce().

Un autre point d’attention est le fait que computeVar() a deux arguments proches mais différents. Le troisième argument –drawsNb– doit être égale à 100 sur chaque noeud mais le quatrième –totalDrawsNb– doit être égal à 1000. totalDrawsNb avec varPrecision permet de définir le nombre de valeurs à retourner: nbValInPercentile = (totalDrawsNb * (1 - varPrecision));. Comme décrit dans mon précédent article, cela permet à chaque noeud de ne retourner que les 10 prix les plus faibles issus de leurs 100 tirages. Cette optimisation de performance réduit le nombre de valeurs envoyées à travers le réseau tout en garantissant que ces 10×100 valeurs envoyées contiendront toujours les 10 plus faibles valeurs des 1000 tirages.

– La fonction reduce() regroupe tous les résultats.

@SuppressWarnings("unchecked")
@Override
public SortedSet<double> reduce(List<gridjobresult> results) throws GridException {
		int nbValInPercentile = 0;
		final SortedSet<double> smallerPrices = new TreeSet</double><double>();

		// Reduce results by summing them.
		for (GridJobResult res : results) {
			SortedSet</double><double> taskResult = (SortedSet</double><double>)res.getData();
			log.debug("Add {} results", taskResult.size());
			smallerPrices.addAll(taskResult);
			Configuration conf = (Configuration)session.getAttribute(SESSION_CONF_KEY);
			if(conf.isCombineAllowed()) {
			//Hypothesis is made that the map task already send only the required values for percentile
				if(nbValInPercentile == 0) {
					nbValInPercentile = taskResult.size();
				}
				else if(nbValInPercentile != taskResult.size()) {
					throw new IllegalStateException("All set should have the same number of elements. Expected " + nbValInPercentile + ", got " + taskResult.size());
				}
				
				while (smallerPrices.size() > nbValInPercentile) {
					smallerPrices.remove(smallerPrices.last());
				}
			}//End if
		}			
		log.debug("Reduce {} jobs and produce {} results", results.size(), smallerPrices.size());
		return smallerPrices;	        	        
}
</double></gridjobresult></double>

Dans notre cas, la tâche reduce identifie les 10 plus petites valeurs dans les 100×10 valeurs les plus faibles envoyées par chaque noeud.

Revenons désormais sur la méthode parametersValueGenerator.intializeRandomGenerator(). Ce paragraphe est un peu technique et peu être ignoré en première lecture. Pour connaître au plus vite comme instancier la grille, vous pouvez aller au paragraphe correspondant.

Problématique de génération de nombres aléatoires

Lorsque vous augmentez le nombre de tirages, la précision devrait augmenter. Durant mes tests, avec un unique thread, la valeur convergeait normalement vers 0.19. Mais avec plusieurs noeuds, la convergence était plus lente et dans certains cas j’obtenais des valeurs surprenantes même avec un fort nombre de tirages. Finalement, après avoir placé des logs à divers endroits du code, j’ai obtenu les lignes suivantes :

14:21:53.635 [gridgain-#17%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 125.53733975440497;122.81494382362781;91.27631322995818;108.45186092917872;126.38622145310778;
14:21:53.635 [gridgain-#4%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 130.49461547615613;113.37411633772678;84.84851329812945;97.68916453583668;126.54608432355151;
14:21:53.635 [gridgain-#6%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 130.49461547615613;113.37411633772678;84.84851329812945;97.68916453583668;126.54608432355151;
14:21:53.636 [gridgain-#7%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 130.49461547615613;113.37411633772678;84.84851329812945;97.68916453583668;126.54608432355151;

La liste SO contient les valeurs générées pour le prix de mon sous-jacent. Vous noterez que différentes séries sont totalement identiques. L’union de ces différents ensembles n’est certainement pas une distribution normale comme nous le souhaitions. Du fait de l’algorithme utilisé et de cette singularité, différents noeuds retournaient des résultats identiques. Comme vous l’avez peut-être noté, ces résultats sont agrégés par une méthode SortedSet.addAll(). Consultez au besoin la JavaDoc de cette méthode et vous découvrirez que les éléments ne sont insérés que s’ils ne sont pas déjà présents. Donc, pour résumé, N noeuds avec des valeurs totalement identiques sont considérés comme équivalents à un unique noeud. Aucun doute, les résultats obtenus seront erronés. Mais pourquoi?

Pour le comprendre, nous devons comprendre comment commons math RandomDataImpl fonctionne. La JavaDoc nous donne un indice en indiquant qu’en réinitialisant avec la même valeur on obtient la même séquence aléatoire « so that reseeding with a specific value always results in the same subsequent random sequence ». L’implémentation par défaut, que j’ai utilisée, se base sur la classe Random du JDK. L’algorithme utilisé pour générer des valeurs aléatoires et appelé générateur congruentiel linéaire. Pour synthétiser, il s’agit d’une série mathématique $X_{n 1}=(a.X_{n} + c) \% n$ qui, étant donné une valeur initiale $X_{0}$ appelée graine (seed en anglais), construit une liste de nombres qui ont les propriétés d’une distribution normale. Pour donner une image concrète, cela peut être comparé – bien que les propriétés mathématiques sous-jacentes soient très différentes- à L’ensemble de Mandelbrot pour lequel une formule mathématique peut conduire à quelque chose d’aléatoire, de chaotique.

(Source http://en.wikipedia.org/wiki/Mandelbrot_set#Image_gallery_of_a_zoom_sequence)
Cependant, si vous donnez deux fois la même racine à un générateur, il générera la même liste de valeurs. Utiliser la fonction initializeRandomGenerator suivante dans le code

public void intializeRandomGenerator() {
	randomGenerator.reSeed(2923455335542820656L);
}

conduit à ce type de log: toutes les listes de valeurs générées sont identiques.

14:51:28.781 [gridgain-#10%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;
14:51:28.781 [gridgain-#14%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;
14:51:28.781 [gridgain-#12%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;
14:51:28.781 [gridgain-#13%null%] DEBUG c.o.r.v.gridgain.VarComputerManager - SO list 105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;105.61456308819832;

Par défaut, la racine utilisée par la classe Random est le timestamp courant. J’ai cherché à utiliser un mélange de l’identifiant de thread et de System.nanoTime() mais j’ai encore obtenu des comportements étranges. La façon la plus simple de corriger le problème a été de définir cet initializedRandomGenerator.

private SecureRandom securedRandomGenerator = new SecureRandom();
//...
public void intializeRandomGenerator() {
		randomGenerator.reSeed(securedRandomGenerator.nextLong());
	}

SecureRandom est conçu pour des scénarios de cryptographie et utilise une implémentation différente qui garantit une distribution plus aléatoire (notamment sans régénérer deux fois les mêmes valeurs). Mais en contrepartie, elle est plus lente que la classe standard Random. Ainsi, cette méthode initializeRandomGenerator() est appelée une fois seulement au début de la méthode computeVar() pour fournir la racine, ensuite la classe standard est utilisée. Cette solution résout le problème. Cependant, il s’agit uniquement d’une première réponse pour les besoins de cet article. Une analyse plus poussée, la comparaison avec d’autres générateurs de nombres aléatoires (cf. cet article Wikipedia) serait nécessaire pour un cas d’utilisation réel.

Instanciation de la grille

Finalement, après avoir résolu ce problème, nous pouvons nous intéresser au code VarComputerGridRunner qui instancie le noeud master de la grille.

public final class VarComputerGridRunner {
	//...
	/**
	 * Method create VarEngine and calculate VAR for a single Call
	 * 
	 * @param args
	 *            Command line arguments, none required but if provided first
	 *            one should point to the Spring XML configuration file. See
	 *            <tt>"examples/config/"</tt> for configuration file examples.
	 * @throws Exception
	 *             If example execution failed.
	 */
	public static void main(String[] args) throws Exception {
		conf = Configuration.parseConfiguration(args);
		// Start GridGain instance with default configuration.
		if (conf.getGridGAinSpringConfigPath() != null) {
			GridFactory.start(conf.getGridGAinSpringConfigPath());
		} else {
			GridFactory.start();
		}

		resultsLog.info("VAR\tDrawsNb\tComputeTime (ms.)\tParameters");
		
		try {
			//Compute multiple times, until 2 min. or 10^6 (32 bits heap size limit)
			double maxTime = System.currentTimeMillis() + 2*60*1000;
			for (int i = 1000; i < 100000001 && System.currentTimeMillis() < maxTime ; i *= 10) {
				VarComputerManager varComputer = new VarComputerManager();
				OptionPricer.Parameters params = new OptionPricer.Parameters(
						252, 120.0, 120.0, 0.05, 0.2);
				OptionPricer op = new OptionPricer(params);
				ParametersValueGenerator valueGenerator = new ParametersValueGenerator(
						0.15, 0, 0);
				printTimeToCompute(varComputer, valueGenerator, op, i);
				valueGenerator = new ParametersValueGenerator(0.15, 0.20, 0);
				printTimeToCompute(varComputer, valueGenerator, op, i);
				valueGenerator = new ParametersValueGenerator(0.15, 0.20, 0.05);
				printTimeToCompute(varComputer, valueGenerator, op, i);
			}

		}
		// We specifically don't do any error handling here to
		// simplify the example. Real application may want to
		// add error handling and application specific recovery.
		finally {
			// Stop GridGain instance.
			GridFactory.stop(true);
			log.info("Test finished");
		}
	}

	private static void printTimeToCompute(
			final VarComputerManager varComputer,
			final ParametersValueGenerator parametersValueGenerator,
			final OptionPricer optionPricer, final int drawsNb)
			throws MathException {
		double startTime = System.nanoTime();
		parametersValueGenerator.initializeReference(optionPricer
				.getInitialParameters());

		final double varPrecision = 0.99;
		SortedSet<Double> smallerPrices = varComputer.computeVar(optionPricer,
				parametersValueGenerator, drawsNb, drawsNb, varPrecision, conf);
		//If combine is not allowed in Map/Reduce, all results are brought back
		//Extraction of percentile should be done here
		if(!conf.isCombineAllowed()) {
			varComputer.extractPercentile(smallerPrices, drawsNb, varPrecision);
		}
		double var = smallerPrices.last();
		double endTime = System.nanoTime();
		resultsLog.info("{}", String.format("%f\t%d\t%f\t%s", var,
				drawsNb, (endTime - startTime) / 1000000, parametersValueGenerator.toString()));
		log.info("{}", String.format("%f (computed with %d draws in %f ms.) with generators %s", var,
				drawsNb, (endTime - startTime) / 1000000, parametersValueGenerator.toString()));
	}
}

Lorsque GridFactory.start() est appelé, le middleware GridGain démarre le noeud master qui recherche les autres noeuds et construit la grille. La configuration par défaut de GridGain utilise le multicast pour la découverte des noeuds mais d’autres implémentations sont possibles. Ensuite, lorsque varComputer.computeVar() est appelé, le middleware GridGain appelle de façon transparente la classe VarComputerGridTask à travers de l’AOP (programmation orientée aspect). La méthode split découpe le calcul, produit plusieurs appels à VarComputer.computeVar() avec différents paramètres. Le middleware envoie le code et les données aux différents noeuds. Chaque noeud exécute alors la fonction VarComputer.computeVar(). Puis les résultats sont renvoyés au noeud maître. Celui-ci exécute la fonction reduce() avant de rendre la main à la fonction main(). Dans mon exemple, la fonction main() est reponsable de :

Identifier le percentile si cela n’a pas été fait par le map/reduce
Identifier la valeur de la VAR : la plus haute valeur du percentile calculé

Quelques mesures de performance

Avec une telle implémentation, il est désormais possible de réaliser quelques mesures de performance. J’ai utilisé mon portable, un DEL Latitude E6510 avec un code i7 quad core avec 4 GB de RAM.

Le premier graphique présente le temps de calcul pour un thread. Il montre que le temps de calcul est proportionnel au nombre de tirages (les deux axes ont des échelles logarithmiques) pour plus de 100000 tirages. Le premier tir avec 1000 tirages est un cas particulier avec un énorme surcoût. Cela est dû au mécanisme de découverte des noeuds. De plus on peut noter que générer 1, 2 ou 3 paramètres n’a pas d’impact visible dans ce cas.

Le second graphique, pour 4 threads simultanés, confirme les points précédents.

Le gain entre 1 et 4 threads n’est pas visible sur ce type de graphiques. Le prochain graphique met en évidence le gain dans deux cas :
$\frac{Compute time N jobs}{Compute time 1 job}$

4 tâches (jobs) lancées dans 4 threads au sein d’un seul processus (Configuration.isThereSeveralSlotsByProcess=true, VarComputerManager est lancé seul)
4 tâches lancées dans 4 processus (Configuration.isThereSeveralSlotsByProcess=false, VarComputerManager est lancé avec 3 autres processus GridGain.bat)

Tout d’abord, et comme on pouvait l’attendre, le gain est meilleur pour un grand nombre de tirages. Ainsi,la distribution ne devrait être utilisée que pour plus de 10000 tirages dans notre exemple.
Ensuite, nous pouvons noter que le gain est inférieur à 4 : entre 2 et 2,5. Les raisons principales que l’on peut évoquer sont :

Mon Core i7 a deux coeurs principaux et la technologie Hyper-Threading autorise 4 threads. Ainsi, même si le système d’exploitation voit 4 processeurs, 2 seulement sont physiquement disponibles. Le gain de performance rendu possible par l’Hyper-Threading est d’après cet article de Wikipedia 30% et non 100% comme nous l’avions supposé en espérant un gain total de 4. Cet article d’Intel (en anglais) va beaucoup plus dans les détails et cite quant à lui des gains jusqu’à 70%. Ainsi 4 threads devraient conduire un gain maximum entre 2×1.3=2.6 et 2×1.7=3.4. L’article d’Intel précise que des applications tournées vers du calcul intensif ont plus de risques de conduire à de plus petits gains avec l’Hyper-Threading. C’est globalement le cas de cette application qui effectue un grand nombre de calculs et utilise peu de mémoire partagée du fait de l’optimisation « combine ». Cette optimisation correspond à envoyer uniquement 1% du nombre total de tirages à la fonction reduce().
Seule la tâche map peut être distribuée. Pour cet algorithme particulier, la fonction reduce est consommatrice en temps
La distribution induit un coût

Enfin, nous constatons que le gain est légèrement plus fort pour les threads avec peu de tirages et plus important pour les processus avec beaucoup de tirages. Comme cette différence n’est pas très marquée, cela peut être dû à une imprécision de la mesure. Si je devais formuler une hypothèse, j’envisagerais que la JVM et l’OS ancrent mieux les processus que les threads sur un coeur permettant une meilleure utilisation du cache.

Pour la suite j’ai utilisé le portable d’un collègue : un Dell Latitude 34800 avec un Core 2 Duo et 2GB de RAM. Les deux machines étaient connectées via notre réseau interne partagé. Il s’agit d’un réseau filaire basé sur un switch. Du fait de la différence du nombre de coeurs, la répartition homogène du nombre de tâches (jobs) par GridGain conduit à 3 tâches sur chaque PC. Ce comportement peut être changé de façon programmatique dans cet exemple. Pour contourner ce problème, j’ai lancé 6 tâches (jobs) dans 6 processus. (Configuration.isThereSeveralSlotsByProcess=false).

VarComputerManager est lancé avec 3 autres processus GridGain.bat processus sur mon portable
2 processus GridGain.bat sont lancés sur l’autre portable

L’impact du premier tir (1000 tirages) est beaucoup plus important. Le temps de calcul est quasi linéaire pour un fort nombre de tirages. Sur ce graphique, les échelles logarithmiques nous empêchent de voir ces gains. Le graphique ci-dessous les rend visibles :

Le Core 2 duo a deux coeurs physiques donc le gain maximum que l’on peut espérer est de 2×1.3+2=4.6. Le gain maximum mesuré se situe aux alentours de 3. Plusieurs raisons peuvent être évoquées:

La distribution sur des machines différentes se fait avec un surcoût plus important que la distribution au sein d’une même machine
La fonction reduce(), qui n’est pas distribuée, est plus chargée car elle doit traiter plus de résultats construits simultanément par les 6 tâches map
Le gain de référence a été mesuré avec 1 thread sur un core i7. Les comparaisons d’Intel effectués avec deux versions différentes du Core i7 et du Core 2 Duo font état de ratio de performance d’environ 1.5. Aussi, une estimation grossière du gain maximum possible dans cette configuration serait de 2×1.3+2×0.66=3.9 ce qui est finalement assez proche de 3.

En somme, la distribution induit un surcoût. Elle doit être limitée aux calculs avec un fort nombre de tirages. Les optimisations faites sur ce code ont permis des gains de performance tout à fait corrects.

Tous les mesures précédentes ont été réalisées avec l’optimisation « combine » activée (consistant à envoyer moins de données à la fonction reduce()). Afin de vous faire sentir le bénéfice de cette optimisation, prenez le temps de regarder les 2 graphiques suivants. Ils montrent le gain relatif à l’optimisation « combine » sur 4 threads sur un PC local et sur 6 processus sur 2 PCs.

Dans les deux cas, utiliser l’option « combine » – i.e. envoyer seulement 10% des données à la fonction reduce() – apporte un gain de X7. Cette valeur est identique dans le deux cas, nous pouvons donc en déduire que l’économie de bande passante réseau n’est pas le point le plus important ici. Les mesures sur mon portable sans option « combine » montre un débit réseau maximal de 11.6 MB/s sur ma carte réseau à 100MB/s. La file d’attente sur cette carte réseau reste à 0 confirmant ce point. Ainsi, nous pouvons en conclure que ce gain a été permis en réduisant la charge de la fonction reduce().

Conclusion

Les middlewares de grille comme GridGain fournissent des fonctionnalités très utiles pour implémenter un calcul de VAR en utilisant la méthode de Monte-Carlo grâce au pattern Map/Reduce. Cependant, malgré leurs efforts, le code doit être conçu en conservant la logique distribuée à l’esprit. Les générateurs de nombres aléatoires doivent être utilisés avec précaution pour rester valides. Partager une configuration ou d’autres informations nécessite d’adapter à la fois le code technique et l’interface du code métier. De plus, des optimisations sont nécessaires pour obtenir de bonnes performances. Diminuer le travail de la fonction reduce() conduit dans mon cas à un gain de x7. Et pour finir, la distribution a un coût. Ce coût est plus important entre machines qu’entre processus locaux. Le surcoût lié à la distribution est rentable pour un grand nombre de tirages. Ainsi la granularité de la distribution doit être optimisée. Cet exemple simplifié fournit un gain relativement bon grâce à la simplification choisie et aux optimisations qui ont pu en découler. Dans un scénario réel, les résultats intermédiaires sont fréquemment stockés afin d’être analysés ultérieurement. Implémenter un tel cas d’utilisation sera le prochain challenge et l’objet de mon prochain article.

2 commentaires sur “Utiliser GridGain pour le calcul de la Value At Risk”

karim aziz
11/05/2011 à 23:18

Bonjour J'ai lu votre article vous êtes vraiment un génie dans le titre (le lien à votre article) vous citez le porcentil(sur exel= porcentil =(tableau, k) k=10% ou 20% etc) je dois calculer en java (mipav) le porcentil des voxel ou pixel) est ce que vous avez un plogin en java s'il vous plait je vous remercie d'avance.

Marc Bojoly
11/05/2011 à 23:38

Bonjour, Pour le calcul d'un percentile en Java je regarderais du côté de Commons Math http://commons.apache.org/math/userguide/stat.html#a1.2_Descriptive_statistics

OCTO respecte votre vie privée

À propos des cookies