Tristan Launay

Anne Philippe

Sophie Lamarche

Construction

Construction of an informative hierarchical prior for a small sample with the help of historical data and application to electricity load forecasting

Keywords: informative prior, hierarchical prior, mcmc algorithms, short dataset, electricity load forecasting

à la diffusion de documents scientifiques de niveau recherche, publiés ou non, émanant des établissements d'enseignement et de recherche français ou étrangers, des laboratoires publics ou privés.

Introduction

We consider the problem of building a prior distribution when only a small number of observations is available. In this context the choice of the prior becomes crucial, having a strong effect on the inference. However it can also be a good tool to bring information so as to improve the estimations in complicated problems (e.g. model in high dimension with respect to the number of observations). Using the ideas of hierarchical modelling, we propose a parametric prior model, that requires the knowledge of another long historical sample y H upon which the model performs equally well in both estimation and prediction. We will assume that this historical sample shares some similarities with our short sample of interest y.

The study is motivated by the modelling and forecasting of electricity load which is a problem well-known within both the academic and the applied statistics community (see e.g. [START_REF] Bunn | Comparative Models For Electrical Load Forecasting[END_REF]. We are interested in the development of a methodology to improve the estimation and the predictions of a parametric multi-equation model (similar to the one presented in [START_REF] Bruhns | A Non-Linear Regression Model for Mid-Term Load Forecasting and Improvements in Seasonnality[END_REF]) over a short dataset. The limited size of the dataset coupled with the high dimensionality of the model leads to a typical overfitting situation when using the maximum likelihood approach: the fitted values are relatively close to the observations while the errors in prediction are an order of magnitude larger or more (due to the very periodic nature of the regressors involved, the model typically requires 4 or 5 years of data to provide satisfactory predictions). This overfitting behaviour can be somewhat alleviated by the use of a Bayesian estimation relying on an informative prior distribution, but the very fact that the data available is limited makes the posterior distribution all the more sensitive to the choice of that prior. Although electricity load curves may largely differ from one population to another, they may also share some common features. The latter case is expected to happen when the global population studied is an aggregation of non-homogeneous subpopulations for which the estimations are made harder due to the relative lack of data.

Different strategies can be found in the literature of hierarchical model to choose or to estimate the parameters of the prior distribution. See [START_REF] Gelman | Bayesian Data Analysis[END_REF]; [START_REF] Gelman | Data Analysis Using Regression and Multilevel/Hierarchical Models[END_REF]; [START_REF] Congdon | Applied Bayesian Hierarchical Methods[END_REF] for a general review on the subject of hierarchical models. These parameters can be assumed unknown and a new level is added to the hierarchical structure. This commonly used approach increases the dimension of the model, which is generally not recommended when dealing with small samples.

In [START_REF] Gelman | Bayesian Data Analysis[END_REF], Chapter 5, some methods using historical samples are presented and illustrated on several examples. Typically, both samples are supposed to be coming from a common model. Assuming that the prior belongs in a parametric family, the hyperparameters of the prior are chosen as the estimated values on the historical sample. This strategy is not robust to small evolutions of the parameters between the historical and the short samples, leading to a bias in the estimation as illustrated in Section 4. Note that it is also possible to estimate the prior on the sample itself, the approaches are often called empirical Bayesian method. This is not a full Bayesian analysis but the approach can be viewed as an approximation to hierarchical models (see e.g. [START_REF] Efron | Large-Scale Inference: Empirical Bayes Methods for Estimation, Testing, and Prediction[END_REF] for details).

Instead of estimating the parameter from historical sample, one can define a joint model on the samples (y H , y) where both samples come from the same parametric model but have different parameters (denoted respectively θ H and θ). The prior distribution on (θ H , θ) conveys the link and the closeness between both parameters with the help of a dependence structure. [START_REF] Bouveyron | Adaptive mixtures of regressions: Improving predictive inference when population has changed[END_REF] propose this method to estimate regression mixtures when the population has changed between the training and the predictive stages. They choose a conditional distribution of θ given θ H which is centered in θ H . This approach requires knowledge of the joint distribution of samples (y H , y), which are assumed independent in [START_REF] Bouveyron | Adaptive mixtures of regressions: Improving predictive inference when population has changed[END_REF]. This leads to a transfer learning in both directions and attributes a higher importance to the largest historical sample.

In this paper, we assume that the long and the short datasets are somehow similar in a nonobvious way. That the similarity between the parameters underlying the two datasets (we will assume they are indeed coming from the model we consider) cannot be easily guessed prevents us from trying to model the datasets simultaneously because it would require a rather precise knowledge of the link between the two. We propose a general way of building an informative hierarchical prior for the short dataset from the long historical one that goes as follows:

1. we first estimate the posterior distribution on the long dataset using a non-informative prior, arguing that the design of an informative prior for this dataset is not necessary, since the data available are enough to estimate and predict the model in this case ;

2. we extract key pieces of information from this estimation to design an informative prior for the short dataset which takes into account the prior information that the datasets are somehow similar, via the introduction of hyperparameters designed to model and estimate this similarity.

The paper is organised as follows. In Section 2 we focus on the general methodology and propose a informative Gaussian prior model. From the two first moments, we show the transfer of information between the historical sample and the prior model. In Section 3 we illustrate the performance of our methodology on simulated dataset from three standard models : the auto-regression model, the polynomial regression and a Poisson model with random effects. In Section 4, we introduce the model used for our main application of electricity load forecasting. We discuss the results obtained on French electricity datasets and show that our new methodology is competitive compared to three alternative standard methods. We also study the effect of the small dataset's sample size upon the predictive quality of the model and show that the informative prior provides reasonable forecasts even when the lack of data is severe. The ad hoc MCMC algorithms we developed to approximate numerically the various estimators are push backed into the Appendix (see Section 5) so as not to obfuscate the main point of the paper by technical details.

Methodology

General principle

We denote y the short sample over which we would like to estimate a parametric model { f θ , θ ∈ Θ} and to predict future observations. We denote y H a historical sample coming from the same parametric model. We assume that this is a long dataset known to share some common features with y. Let θ H and θ be the corresponding parameters. The parameter of interest is θ .

Let π H and π be the prior distributions used to estimate θ H and θ respectively. Let π H (•|y H ) and π(•|y) be the associated posterior distributions. We will denote θ H and θ the corresponding Bayes estimates (i.e. the posterior means).

We propose a method to the improve parameter the estimations and model predictions for the dataset y with the help of the dataset y H . Observe first that the choice of π H is not crucial as long as it remains non-informative enough because, θ H can be correctly estimated from the data alone on y H . The naive pick π = π H (•|y H ) is not a viable solution as soon as the parameters of θ H and θ differ since the variance of the posterior distribution π H (•|y H ) is too small: in that case the few observations of y will not be able to make up for the difference between the posterior mean θ H and the true value of the parameter θ. Assuming that the parameters θ H and θ are identical is too restrictive in practise. To allow for more flexibility we introduce hyperparameters to account for the similarity between the datasets. We now describe the informative hierarchical prior we designed.

Assume that the prior distribution π is to be chosen within the parametric family

F = {π λ ; λ ∈ Λ}.
Since selecting π ∈ F is equivalent to picking λ ∈ Λ, and since we would like π to retain some key-features of π H (•|y H ), we want to pick λ using some of the information contained inside the posterior distribution π H (•|y H ). We assume that there exists an operator T : F → Λ, such that

T[π λ ] = λ,
and choose λ proportional to T[π H (•|y H )], in the sense that λ = KT[π H (•|y H )],
where K : Λ → Λ itself is an unknown linear operator that we assume diagonal for ease of use.

The operator K can be interpreted as a similarity operator between θ H and θ, and its diagonal components as similarity coefficients measuring how close the two datasets y H and y really are when looked at through the model. The diagonal components of K are hyperparameters of the prior we designed, and we give them a vague hierarchical prior distribution centred around q, the prior on q being vague and centred around 1.

The hyperparameter q may also be regarded as a more global similarity coefficient, since it represents the mean of all the similarity coefficients. The prior mean of q is forced to 1 to reflect the prior knowledge that the datasets are somehow similar. The variance of the prior distribution of q could in theory be reduced, going from a vague prior to a more informative structure, depending on the confidence we have over the similarity between the datasets. We chose not to however, so as to keep the procedure we describe from requiring any delicate subjective adjustments.

Remark 1. As we suppose that y H is a large sample, the posterior distribution π H (•|y H ) can be approximated by Gaussian distribution under suitable regularity conditions. Therefore we can replace

λ = KT[π H (•|y H )], by λ = KT[N ( θML H,n , I -1 H,n )],
where θML H,n is the maximum likelihood estimate and where I H,n is an estimate of the expected Fisher Information matrix calculated on the sample y H (Ghosh et al., 2006, chap 4).

We present now two frequent situations where the above procedure can be written in a simpler way.

Example 2 (Method of Moments). We assume that the elements of F can be identified via their m first moments: the operator T can then be reduced to a function F of the m first moments operators, i.e.

λ = T[π λ ] = F(E(θ), . . . , E(θ m )). The expression of λ then becomes λ = KF(E(θ|y H ), . . . , E(θ m |y H )).
Note that, if the prior requires the specification of at least the two first moments, even though the priors from the upper layers of the model are vague, the correlation matrix estimated on the dataset y H remains untouched and is directly plugged into in the informative prior if we consider centred moments for orders greater than 1.

Example 3 (Conjugacy). We consider the case where F is the family of priors conjugated for the model. If the prior π H belongs to F then the associated posterior distribution π H (•|y H ) does too and we denote λ H (y H ) the parameter that corresponds to it. The expression of λ then reduces to λ = Kλ H (y H ).

Gaussian strategy

Under suitable conditions of regularity (see [START_REF] Ghosh | An Introduction to Bayesian Analysis, Theory and Methods[END_REF], [START_REF] Launay | Consistency of the posterior distribution and MLE for piecewise linear regression[END_REF]), the posterior distribution tends to look like a normal distribution asymptotically. This suggests the use of a normal distribution as the prior distribution π. Let us denote µ H and Σ H the posterior mean and posterior variance of the parameter θ H .

As described in the general methodology, we introduce new parameters to model the similarity between the two datasets, (k,

l) ∈ R d × R and (q, r) ∈ R × R * + . The informative hierarchical prior is of the form θ|k, l ∼ N (K µ H , l -1 Σ H ) (1) k|q, r ∼ N (q(1, . . . , 1) ′ , r -1 I d ) (2)
where K = diag(k). Again, the coordinates of the vector k can be interpreted as similarity coefficients between θ H and θ and the strictly positive scalar l can be seen as a way to alternatively weaken or strengthen the covariance matrix as needed. Hyperparameters q and r are more general indicators of how close θ H and θ are, q corresponding to the mean of the coordinates of k and r being their inverse-variance. For a complete definition of the prior distribution, we of course need to specify the prior on l, q, r. Based on conjugacy considerations, we use:

l ∼ G(a l , b l ), q ∼ N (1, σ 2 q ), r ∼ G(a r , b r ), (3)
where a l , b l , a r , b r and σ 2 q are fixed positive real numbers such that the prior distribution on l, q and r are vague. These prior distributions are chosen because of their conjugacy properties (that we took advantage of to design the MCMC algorithms for our applications). The vagueness requirement that we impose on these priors is motivated by the fact that we want to keep as general a framework as possible without having to tweak each and every prior coefficient for different applications.

In the end, the informative prior is built as follows:

π(θ, k, l, q, r) ∝ π(θ|k, l)π(k|q, r)π(l)π(q)π(r) (4)
where the densities appearing on the right hand side are specified in (1), ( 2) and (3).

From the two first moments of the prior distribution π, we evaluate the prior information brought by the historical sample y H . Firstly we show that the prior distribution π is centered on the estimation of the parameter on y H , which should be close to θ H according to the asymptotic properties of Bayes estimates.

Property 4. The prior distribution admits a moment of order 1 equal to µ H . Proof. Indeed, we have

E(θ) = E(E(θ|K, l -1 )) = E(K µ H ) = E(K) µ H = E(E(K|q, r)) µ H = E(q)I d µ H = µ H .
Secondly, from the correlation matrix of π, we show that we transfer informations from the historical sample to the dependence between the components of parameter θ Property 5. Assume that a l > 1 and a r > 1. Under suitable conditions of regularity,

• the correlation of the prior distribution is equivalent to the correlation of the posterior distribution π H (•|y H ) as N H goes to infinity.

• for a fixed N H the correlation of the prior distribution is also equivalent to the correlation of the posterior distribution π H (•|y H ) as σ 2 q and E(r -1 ) go to 0.

Proof. First notice that the condition on the coefficients a l > 1 and a r > 1 ensures that the first moment of both l -1 and r -1 are finite. This implies that the variance of the prior distribution is finite Recall now that if (X, Y) are two random variables in L 2 we have the following equality var(X) = var(E(X|Y)) + E(var(X|Y)).

We then get var(θ) = var(E(θ|K, l -1 )) + E(var(θ|K, l -1 ))

= var(K µ H ) + E(l -1 ) ΣH = diag( µ H ) var(k)diag( µ H ) + b l a l -1 ΣH .
Using the same equality, we can evaluate the variance of k as var(k) = var(E(k|q, r -1 )) + E(var(k|q, r -1 ))

= var(q    1 . . . 1   ) + E(r -1 )I d = σ 2 q    1 • • • 1 . . . . . . . . . 1 • • • 1    + b r a r -1 I d .
We obtain

var(θ) = diag( µ H )   σ 2 q    1 • • • 1 . . . . . . . . . 1 • • • 1    + b r a r -1 I d    diag( µ H ) + b l a l -1 ΣH = b l a l -1 Σ H + b r a r -1 diag( µ H ) 2 + σ 2 q µ Ht µ H .
and the element of the correlation matrix in i th row and j th column

(i = j) is thus b l a l -1 Σ H (i, j) + σ 2 q µ H (i) µ H (j) b l a l -1 Σ H (i, i) + ( b r a r -1 + σ 2 q ) µ H (i) 2 1/2 b l a l -1 Σ H (j, j) + ( b r a r -1 + σ 2 q ) µ H (j) 2 1/2 .
Under suitable conditions of regularity, µ H converges to θ H and 1

N H

ΣH converges to the Fisher information I(θ H ) as N H tends to infinity, which ends the proof of the first assertion. The second assertion is trivial when you recall that E(r -1 ) = b r a r -1 .

Basic examples

In this Section, we describe how to apply our methodology to three different and common examples. We first deal with both a regression model and an auto-regressive model, before addressing a hierarchical Poisson model.

Regression model

We will consider a Gaussian regression model i.e. for each t = 1, . . . , n

y t = f t (θ) + ǫ t (5)
where ǫ 1 , . . . , ǫ N are assumed independent and identically distributed with common distribution N (0, σ 2 ).

We will use the notation f (θ) = ( f 1 (θ), . . . , f N (θ)) for short, where θ designate the parameters of interest. With these notations, since η = (θ, σ 2 ), the likelihood of the model considered is

L(y|η) = √ 2πσ -N exp - 1 2 y -f (θ) 2 . ( 6
)
For σ 2 we use a non-informative prior (we chose π(σ 2 ) ∝ σ -2 ) because we do not want to make any kind of assumptions about the noise around both datasets. This prior is noninformative in the sense that it matches Jeffreys' prior distribution on σ 2 for a Gaussian linear regression. For the parameter θ we use the informative hierarchical prior as described in Section 2.2

On the sample y H we use the following non-informative prior

π H (η H ) = π H (θ H , σ 2 H ) ∝ σ -2 H . ( 7
)
This prior is non-informative in the sense that it matches Jeffreys' prior distribution on σ 2 H for a Gaussian linear regression (see Marin and Robert, 2014, on page 73). It leads to the following posterior distribution

π H (η H |y) ∝ L(y H |η H )π H (η H ) (8) ∝ σ -N H -2 H exp - 1 2 σ -2 H y -f (θ H ) 2 2 .
We consider a polynomial regression model of order p (p ≥ 0) on the interval [-1, 1], meaning that observations from this model can be written for each t = 1, . . . , n

y t = p ∑ i=0 θ i x p t + ǫ n (9)
where ǫ n are i.i.d. with Gaussian distribution N (0, σ 2 ), and where x 1 , . . . , x n are equally spaced over [-1, 1].

Here again, we only show some of the results we obtained in various settings. Our goal here is to show, for a fixed sample size, the robustness of our method to non-exact similarity between the two datasets First situation We simulate N H = 200 observations y H from model ( 9) with p = 4 and coefficients θ H = (2, -1, 3, 1, 2)

and N = 10 observations y with coefficients θ = ρθ H , where ρ will be chosen later, using σ 2 = 1 for the variance of the noise for y H and σ 2 = 4 for y. We then estimate the coefficients θ from the small sample y available using the Gaussian informative prior we presented above. To assess the performance of our method, we compare and aggregate the results over 100 repetitions (sampling new observations for each repetition) with those obtained from a non-informative vague prior and those coming from MLE estimations.

Second situation

The second situation is different from the first in that θ is not proportional to θ H anymore. Here we chose

θ 0 = ρθ H 0 θ 1 = ρθ H 1 θ 2 = θ H 2 θ 3 = θ H 3 θ 4 = θ H 4
and the rest of the setup is identical.

Results. Figure 1 shows the errors of estimation for the four coefficients i.e. the true values minus the estimated values (the posterior mean for Bayesian models) when ρ varies from 1 to 0.5 for both situations and Figure 2 displays the actual estimated polynomials for both situations when ρ = 0.5. Once again the benefits of using the informative prior over the non-informative prior stand out. Whether the similarity between the two parameters θ and θ H is based upon proportionality or not the informative prior helps the estimation of all the coefficients in a very effective way.

As the sample size N decreases from 100 to 10, the MLE and the non-informative Bayesian method tend to perform rather poorly whereas the use of the informative Bayesian prior clearly benefits the estimation of the coefficients in this situation where observations y H and y come from exactly the same model.

Auto-regressive model

We consider an AR(p) model (p > 0) with coefficients θ, and variance of noise σ 2 , i.e. denoting y n its observation at time n, for each n

y n = p ∑ i=1 θ i y n-i + ǫ n ( 10
)
where ǫ n are i.i.d. with Gaussian distribution N (0, σ 2 ). As in the previous example, we choose for θ the informative hierarchical prior described in Section 2.2, and for σ 2 we take the non-informative prior proportional to σ -2 . On the historical sample y H , we again choose the non informative prior (7) (see Marin and Robert, 2014, on page 219).

We tested our method on many different such models (with varying orders, coefficients, and variance of noise) and present below only the results obtained for an auto-regressive model of order 4 (p = 4) and show the effect of the sample size N.

We simulate N H = 400 observations y H from model (10) with coefficients θ H = (1.7, -0.23, -0.833, 0.3528) (the corresponding polynomial having roots -0.7, 0.7, 0.8 and 0.9), and also simulate N = 100, 30, 10 observations y from the same AR(4) model with coefficients θ = θ H , using in each case σ 2 = 1 for the variance of the noise. We then estimate the coefficients θ from the observations y available using the Bayesian prior we presented above. To assess the performance of our method, we compare and aggregate the results over 100 repetitions (sampling new observations for each repetition) with those obtained from a non-informative vague prior and those coming from MLE estimations. Figure 3 shows the errors of estimation for the four coefficients i.e. the true values minus the estimated values (the posterior mean for Bayesian models).

Hierarchical Poisson model

We consider y 1 , . . . , y N independent random variables such that for each j, y j follows a Poisson distribution with parameter λ j . We assign a Gamma distribution for the parameters λ i : When a large sample is available, we can choose a non informative distribution on the parameter θ = (a, b), for example the Jeffreys prior which is given by [START_REF] Yang | A catalog of noninformative priors[END_REF] ). When a goes to zero or to infinity, a(log(Γ(a)) ′′ -1 is equivalent to 1 √ a (see [START_REF] Abramovitz | Handbook of Mathematical Functions with Formulas, Graphs, and Mathematical Tables[END_REF] chapter 6). To simplify the computational aspects, we consider a prior distribution which behaves as the Jeffreys prior when a tends to zero and to infinity. We assume that a and b are independent, the marginal distributions of a and b are Gamma distributions with parameters (1/2, ǫ) and (ǫ, ǫ) respectively. The parameter ǫ is fixed small enough.

y i ∼ P (λ i ) λ i ∼ Γ(a, b)
π J (a, b) = 1 b a(log(Γ(a)) ′′ -1 (see
For this choice of prior, Figure 4 illustrates that the Bayes estimates do not perform well for the small sample.

In the context of Section 2, there is a large sample y H sharing some similarities with the small sample of interest y. Therefore we can construct the prior distribution on (a, b) using the sample y H based on either of the two following strategies 1. we assign on (a, b) the hierarchical prior defined in (4) with θ = (a, b). Figure 4 shows the improvement brought by the informative hierarchical model upon non informative model. On independent replications, we estimate the mean square error on the estimation of

θ i E N = 1 N N ∑ i=1 |λ i -E(λ i |y)| 2 . ( 11
)
The results are given in Figure 5. When the two populations have exactly the same parameters, it is not surprising to obtain better results using fixed parameters. However the simulation shows also that this strategy is not robust and leads to worse results than the non informative prior when the true parameters differ. The hierarchical prior on (a, b) gives some very satisfying results in terms of quality of the estimation and robustness to variation of the parameters between the sample y H and y.

• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 25
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 25
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 25
• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 25
(a H , b H ) = (a, b) = (1, 1)
. The boxplot is constructed from 500 independent replications.

• • • • • • • • • • • • • • • • • • • • • • • • • Info. Fixed ni Info. Fixed ni Info. Fixed ni 0.0 0.5 1.0 1.5 2.0 n=25 • • • • • • • • • • • • • • • • • • • • • Info. Fixed ni Info. Fixed ni Info. Fixed ni 0.0 0.5 1.0 1.5 2.0 n=50 • • • • • • • • • • Info. Fixed ni

Electricity load forecasting

Modelling and forecasting the electricity load on a day-to-day basis has long been a key activity for any company involved in the electricity industry. It is first and foremost needed to supply a fixed voltage at all ends of an electricity grid: the amount of electricity produced has to match the demand very closely at any given time and experts usually make use of short-term forecasts (see [START_REF] Cottet | Bayesian modeling and forecasting of intraday electricity load[END_REF], for example)).

Electricity load usually involves a large predictable component due to its very strong daily, weekly and yearly periodic behaviour and the exogenous variables most commonly used to forecast the electricity load are weather-based (see [START_REF] Al-Zayer | Modelling the Impact of Temperature on Electricity Consumption in the Eastern Province of Saudi Arabia[END_REF][START_REF] Taylor | Using weather ensemble predictions in electricity demand forecasting[END_REF], for an example). For the French electricity load specifically, the importance of the temperature and cloud-cover is underlined in [START_REF] Menage | Nouvelle modelisation de l'influence des conditions climatiques sur la consommation d'energie electrique[END_REF] and [START_REF] Bruhns | A Non-Linear Regression Model for Mid-Term Load Forecasting and Improvements in Seasonnality[END_REF]. [START_REF] Engle | Semiparametric estimates of the relation between weather and electricity[END_REF]; [START_REF] Ramanathan | Short-run forecasts of electricity loads and peaks[END_REF] also present other weather-related models. Some authors work with univariate time series models: [START_REF] Taylor | Short-term electricity demand forecasting using double seasonal exponential smoothing[END_REF] builds a double seasonal exponential smoothing for the British electricity load, and [START_REF] Cugliari | Prévision non paramétrique de processus à valeurs fonctionnelles, application à la consommation d'électricité[END_REF] treats the load curve as a functional-valued autoregressive Hilbertian process to use a non-parametric approach relying on the wavelet transform. Other authors (see [START_REF] Soares | Modeling and forecasting short-term electricity load: a comparison of methods with an application to brazilian data[END_REF][START_REF] Ramanathan | Short-run forecasts of electricity loads and peaks[END_REF][START_REF] Dordonnat | An Hourly Periodic State Space Model for Modelling French National Electricity Load[END_REF], for example) choose multiple-equation models where the various hours of the day share the same equation, while the associated parameters differ.

Bayesian methods have also been used in the past to solve electricity load forecasting problems. Amongst them, [START_REF] Minka | Bayesian linear regression[END_REF] deals with Bayesian linear regression and [START_REF] Harrison | Bayesian forecasting[END_REF] focus on dynamic linear models from the Kalman filter point of view and [START_REF] Smith | Modeling and short-term forecasting of new south wales electricity system load[END_REF] choose a Bayesian semi-parametric regression model. To the best of our knowledge though, no previous attempt has been made in the literature to help improve the forecasting of an intraday electricity load model on a short dataset, based upon another dataset using a hierarchical Bayesian model.

Regression model

In the next few paragraphs, we specify the regression model that we use to model the electricity load.

For each of the 48 instants of the day (each instant lasts 30 minutes, starting from 00:00am), the non-linear regression model that we consider in this paper, first described in [START_REF] Bruhns | A Non-Linear Regression Model for Mid-Term Load Forecasting and Improvements in Seasonnality[END_REF], is made of three components and is usually formulated as follows: for t = 1, . . . , N,

y t = x (1) t x (2) t + x (3) t + ǫ t (12)
x

(1)

t = d 11 ∑ j=1 z cos j cos 2jπ 365.25 t + z sin j sin 2jπ 365.25 t + d 12 ∑ j=1 ω j 1 Ω j (t), x (2)
t = d 2 ∑ j=1 ψ j 1 Ψ j (t), x (3)
t = g(T t -u)1 [T t , +∞[ (u),
where y t is the load of day t and where ǫ 1 , . . . , ǫ N are assumed independent and identically distributed with common distribution N (0, σ 2 ). The x (1) component is meant to account for the average seasonal behaviour of the electricity load, with a truncated Fourier series (whose coefficients are z cos j ∈ R and z sin j ∈ R) and gaps (parameters ω j ∈ R) which represent the average levels of electricity load over predetermined periods given by a partition (Ω j ) j∈{1,...,d 12 } of the calendar. This partition usually specifies holidays, or the period of time when daylight saving time is in effect i.e. major breaks in the electricity consumption behaviour. The left part of Figure 6 shows the typical behaviour of the electricity load over two different periods of time (summer vs. winter). The x (2) component allows for day-to-day adjustments of the seasonal behaviour x (1) through shapes (parameters ψ j ) that depends on the so-called days' types which are given by a second partition (Ψ j ) j∈{1,...,d 2 } of the calendar. This partition usually separates weekdays from weekends, and bank holidays. The differences between two different daytypes are visible on the left part of Figure 6 too. For obvious identifiability reasons, the vector ψ is restricted to the positive quadrant of the • 1 -unit sphere in R d 2 , that we denote

• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • ••• • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • ••• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • •• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • -5 0 5
S d 2 + (0, 1) = ψ ∈ (R + ) d 2 ; ψ 1 = 1 .
The x (3) component represents the non-linear heating effect that links the electricity load to the temperature (see [START_REF] Seber | Nonlinear Regression[END_REF], for a general presentation of non-linear models), with the help of 2 parameters. The heating threshold u ∈ [u, u] corresponds to the temperature above which the heating effect is considered null and is usually estimated to be roughly around 15 • C. The heating effect is supposed to be linear for temperatures below the threshold and null for temperatures above. Note that the threshold is sought within the range of the observed temperatures, i.e. u ∈ [u, u] with min t=1,...,N T t < u < u < max t=1,...,N T t .

The heating gradient γ ∈ R * where R * = R\{0} represents the intensity of the heating effect, i.e. the slope (assumed to be non-zero) of the linear part that can be observed on the right part of Figure 6.

Priors

For the sake of both clarity and completeness, we present the full Bayesian model in use for all our applications in the lines below. Let us denote here again η = (θ, σ 2 ) where θ designates the parameters of interest in the electricity load model (i.e. z j , ω j , ψ j , g, and u).

We opted for the Gaussian strategy (see Section 2.2) for the regression model as was described in Section 3.1, meaning that the non-informative prior distribution π H used with the historical sample y H is

π H (η H ) = π H (θ H , σ 2 H ) ∝ σ -2 H .
The informative prior π used with the small sample y is defined hierarchically and can be summarised as π(θ, σ 2 , k, l, q, r) ∝ π(θ|k, l)π(k|q, r)π(l)π(q)π(r)π(σ 2 ) with:

θ|k, l ∼ N (K µ H , l -1 Σ H ) q ∼ N (1, σ 2 q ) k|q, r ∼ N (q(1, . . . , 1) ′ , r -1 I d ) r ∼ G(a r , b r ) l ∼ G(a l , b l ) π(σ 2 ) ∝ σ -2 ,
where K = diag(k), where µ H and Σ H are the posterior mean and posterior variance of the parameter θ H and where a l , b l , a r , b r and σ 2 q are fixed positive real numbers such that the prior distribution on l, q and r are vague. Remark 6. Using either the non-informative prior distribution or the informative prior distribution with the model (12) leads to a proper posterior distribution. The proofs, as well as the MCMC algorithms we designed to derive numerical estimators, are pushed back into the Appendix (see Section 5).

Datasets

The historical dataset y H needed for the construction of the informative prior corresponds to a specific population in France frequently referred to as "non-metered" because their electricity consumption is not directly observed by EDF but instead derived as the difference between the overall electricity consumption and the consumption of the "metered" population. For this population the data ranged from 07/01/2004 to 07/31/2010.

We illustrate the benefits of choosing our informative prior to predict electricity load on two short datasets : the first dataset y corresponds to the "non-metered" population for ERDF, a wholly owned subsidiary of EDF that manages the public electricity network for 95% of continental France. This population roughly covers the same people that y H does, but not exactly. The population behind the second dataset y ′ corresponds to a subset of the population behind y H and represents around 50% of the total load of y H .

Benchmark against standard methods

For this application, we only consider the days for which no special tariffs are enforced: the so-called EJP ("Effacement jour de pointe" = peak tariff days) are removed from the dataset beforehand to ensure the signal studied is consistent throughout time. Bank holidays (including the day before and the day after to avoid any neighbourhood contamination effects), the summer holiday break (August) and the winter holiday break (late December) are also removed from the dataset for this first application, so that we may benchmark our method against others on a smooth and rather easy-going signal, and not put any one method at a disadvantage due to the signal's specificity. The temperature considered in the model is the average temperature over France for the period of study, and the cooling threshold is chosen to be 16 • C throughout the 48 instants of the day. We benchmark our Bayesian method with informative prior against four alternative methods, comparing their predictions on dataset y in two configurations. Roughly speaking, for our first experiment (Case 1) we estimate the model for y over the period ranging from 12/01/2009 to 06/30/2010 and predict the next 30 days (same as the application presented Section 4.5), while for our second experiment (Case 2) we estimate the model for y over the period ranging from 01/01/10 to 07/31/10 and predict the previous 30 days. We expect Case 1 to be the easy case and Case 2 to be the tough case, the signal being very smooth during summer and not so much during winter. The figures shown in Table 1 summarise the exact lengths of the various datasets for both experiments.

The four alternative methods we benchmark against our own Bayesian informative method, rely on four different techniques: the first one is the Bayesian non-informative method that we presented earlier in Section 2.2 (recall that it is meant to be an equivalent to the maximum likelihood approach), the second involvedsnon-parametric estimation with kernels (see [START_REF] Fan | Non linear Time Series: Nonparametric and Parametric Methods[END_REF], the third is a double exponential smoothing (see [START_REF] Taylor | Short-term electricity demand forecasting using double seasonal exponential smoothing[END_REF] and the fourth and last is an ARIMA model. Note that for the second experiment the data available obviously has to be time-reversed in order to apply the last three alternatives methods, since time-dependence plays an important role for them. The ARIMA model is automatically selected (see [START_REF] Hyndman | Automatic time series forecasting: the forecast package for R[END_REF] as the best model with regard to the AIC criterion.

It is clear from the results exposed in Table 2 that the informative prior outperforms all the alternative methods by a large margin in each case. The overall bad performance of the Bayesian non-informative method is not surprising because at least 3 to 4 years of data are usually required to avoid overfitting, for such a parametric model.

Role of the hyperparameters

For this application, the setup is the same as the one described at the start of Section 4.4. Our aim here is to point out the role of the hyperparameters introduced within the informative prior and show that besides providing better results than alternative methods as was demonstrated in Section 4.4, they also provide a measure of similarity between the short datasets of interest and the dataset used to build the prior. We estimate the model for y and y ′ over the period ranging from 12/01/2009 to 06/30/2010 and predict the next 30 days.

Figure 7 shows the predictive quality of the model for both populations y and y ′ using the informative prior as well as the posterior densities of the similarity coefficients k j at midday (we also checked the other 47 instants but the look of them was nearly identical to the one we present). These densities are much more peaked and also centred closer to 1 for population y than they are for population y ′ . It thus seems to indicate that dataset y is more similar to y H than y ′ is, confirming our prior knowledge that y H and y cover approximately the same population whereas y ′ represents around 50% of y H : this value of 50% is also visible on Figure 7 where we observe two densities centred around 0.5, which correspond to the similarity coefficients between the offsets ω j of y H and these of y ′ .

Figure 8 displays the boxplots for the posterior densities of q and 1/ √ r and seem to corroborate the fact that y is more similar to y H than y ′ is. Recall that q and 1/ √ r respectively act as the mean and standard deviation of the similarity coefficients k j within our informative hierarchical prior. Indeed the estimated mean of q appears to be closer to 1 while its estimated variance is smaller on y than y ′ . The estimated mean and variance of 1/ √ r are also smaller on y than y ′ . The estimated values of q and r provide some information about the similarity between the datasets considered. Notice also that the best predictive performance is obtained when the similarity between the two datasets is strongest.

Role of the sample size

For this application, the setup is almost identical to the one described at the start of Section 4.4 but the temperature used in the model is not the average temperature over France anymore but a transformation of it: it is smoothed using exponential smoothing, which is known to improve the link existing between the two variables temperature and electricity load (see [START_REF] Bruhns | A Non-Linear Regression Model for Mid-Term Load Forecasting and Improvements in Seasonnality[END_REF] for more information about this). The cooling threshold is fixed at 18 • C throughout the 48 instants of the day, and this time, the summer holiday break is not removed from the dataset (but the winter holiday break and the bank holidays still are), so that the model could benefit from (and be tested on) the August months in general. Note that for y H we use the same dataset as we used for our two first applications.

We put the focus on the length of the estimation period on y while keeping the same prediction window. We successively choose the periods ranging from 01/01/2010, 03/01/2010, 05/01/2010, 07/01/2010 to 12/31/2010, reducing the estimation period on y from 12 months to only 6 months, removing 2 months at a time. We then predict the next 6 upcoming months i.e. from 01/01/2011 to 06/30/2011. The diagram in Figure 9 describes the 4 scenarios considered.

The non-informative prior leads to a better fit than the informative prior as can be seen in Table 3. It should not come as a surprise because the non-informative prior was indeed meant to be equivalent to a maximum likelihood approach whose criterion is precisely to minimise the RMSE. As for the quality of the predictions associated with the model for both priors, Table 3 demonstrates that the informative prior beats the non-informative prior in each of the four proposed configurations. The improvement appears to be minimal when 12 months are used but as months are removed from the estimation window, the predictive quality for the non-informative prior drops very quickly, while the predictive quality for the informative prior remains moderate and stable.

Figure 10 shows the average error in prediction for each month. It is important to note that the use of the non-informative prior leads to overfitting the model: the results presented in Table 3 show that as the estimation window goes smaller, the estimation error decreases while the prediction error grows very quickly. A close inspection of the posterior densities of the different parameters of the model reveals that the bias induced by the increasing lack of data is mainly seasonal: this is due to the seasonality coefficients of the model being overfit. Choosing the informative prior over the non-informative prior makes the estimation and prediction of the Figure 8: Boxplots of the posterior densities of q (left) and 1/ √ r (right) (mean and standard deviation of the similarity coefficients k j ) for populations y and y ′ , throughout the 48 instants of the day. model more robust with regards to the lack of data.

• • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • 0 10 20 30 40 2 4 6 8 10 MAPE • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • • B B
The informative prior especially improves the quality of the predictions when the lack of data is severe: it provides reasonable forecasts even in the worst scenario considered here, where only 6 months of estimation are used for 6 months of prediction. In this situation, estimation (from 07/01/2010 to 12/31/2010) and prediction (from 01/01/2011 to 06/30/2011) are performed on non-overlapping areas of the calendar: the informative prior makes up for the unavailable data and prevents the model from overfitting on the second half of the calendar while the non-informative prior does not and consequently leads to heavily biased predictions over the first half of the calendar.

Jan.

Feb. Mar.

Apr. May Jun. Figure 10: Quality of the predictions (MAPE in %) averaged for each month (all the instants of the 30 or so days within each month are used to compute these averages), with an estimation period ranging from 12 to 6 months using the non-informative prior (left) and using the informative prior (right). The ordinate axis is in log-scale. Each shade of grey corresponds to a different scenario (as depicted on Figure 9).

Non-informative prior (NI)

Let (β 0 , u 0 ) ∈ B d β + (0, 1) × [u, u] and denote α * = (α, γ). We write y -f ((α, β 0 , γ, u 0 )) 2 2 = N ∑ t=1 y t -(B t• β 0 + C t )A t• α -(T t -u 0 )1 [T t , +∞[ (u 0 )γ 2 = y -A * (β 0 , u 0 )α * 2 2 ,
and thus obtain the following equivalence, as (β, u) → (β 0 , u 0 ) and α

* 2 → +∞ y -f (θ) -N 2 ∼ y -A * (β 0 , u 0 )α * -N 2 . ( 14
)
The triangular inequality applied to the right hand side of ( 14) gives

y -A * (β 0 , u 0 )α * -N 2 ≤ y 2 -A * (β 0 , u 0 )α * 2 -N . ( 15
)
Since A ′ * (β 0 , u 0 )A * (β 0 , u 0 ) has full rank, by straightforward algebra we get

λ α * 2 2 ≤ A * (β 0 , u 0 )α * 2 2 ,
where λ is the smallest eigenvalue (A * (β 0 , u 0 )) ′ A * (β 0 , u 0 ) and is strictly positive. We can hence find an equivalent of the right hand side of (15) as α * 2 → +∞, which is

y 2 -A * (β 0 , u 0 )α * 2 -N ∼ λ -N/2 α * -N 2 . ( 16
)
Combining ( 14), ( 15) and ( 16) together, we see that the integrability of the left hand side of ( 14) as (β, u) → (β 0 , u 0 ) and α * 2 → +∞ is directly implied by that of α * -N

2 . The latter is immediate for N > d α + 1, as can be seen via a Cartesian to hyperspherical re-parametrisation.

The previous paragraph thus ensures the integrability of yf (θ) -N 2 over sets of the form

{(β, u) ∈ V((β 0 , u 0 )), α * 2 ∈]M(β 0 , u 0 ), +∞[}, ∀(β 0 , u 0 ) ∈ B d β + (0, 1) × [u, u]
where the subset V((b 0 , u 0 )) is an open neighbourhood of (β 0 , u 0 ) and M(β 0 , u 0 ) is a real number depending on (β 0 , u 0 ). By compacity of B d β + (0, 1) × [u, u] there exists a finite union of such V((β i , u i )) that covers B d β + (0, 1) × [u, u]. Denoting M the maximum of M(β i , u i ) over the corresponding finite subset of (β i , u i ), we finally obtain the integrability of y -

f (θ) -N 2 over {(β, u) ∈ B d β + (0, 1), α * ∈]M, +∞[}. The integrability of y -f (θ) -N 2 over {(β, u) ∈ B d β + (0, 1), α * ∈ [0, M]} is trivial, recalling that θ → y -f (
θ) 2 is continuous and does not vanish over this compact for almost every y, meaning its inverse shares these same properties.

Remark 9. The condition "A ′

* A * has full rank" mentioned above is typically verified in our applications for the regressors used in our model. To see this, call "vector of heating degrees" the vector whose coordinates are (T tu)1 [T t , +∞[ (u), then not verifying the aforementioned condition is equivalent to saying that there exists an index i and a threshold u such that the family of vectors formed by the regressors A and the vector of heating degrees is linearly dependent over the subset Ψ i of the calendar".

MCMC algorithms for the estimation of the posterior density

The two MCMC algorithms presented below were developed because direct simulations from the posterior distribution were not possible. The justifications are given after the algorithms themselves. Notice that the full conditional distributions of all the parameters but the threshold u appear to be common distributions in both cases, due to the presence of multiple semi-conjugacy situations. We used a Metropolis-within-Gibbs algorithm (see Robert and Casella, 2009, page 230, for a description) based on Gibbs sampling steps for every parameter but u for which we used a Metropolis-Hasting step based on a Gaussian random walk proposal.

MCMC algorithm for the estimation of the posterior distribution, using the informative prior

In the lines below we give the different steps of the MCMC algorithm we used to (approximately) simulate (η 1 , . . . , η M ) according to the posterior distribution π(η|y) corresponding to the informative prior we presented earlier. Denoting N (µ, Σ, S) the S-truncated Gaussian distribution with mean µ and covariance Σ, the algorithm goes as follow:

Step 1. Initialise η 1 such that π(η 1 |y) = 0

Step 2. For t = 1, . . . , M -1, repeat (i). Simulate σ 2 t+1 cond. to (α t , β t , γ t , u t , k t , l t , q t , r t , y)

σ 2 t+1 ∼ IG N 2 , 1 2 y -f (θ) 2 2
(ii). Simulate r t+1 cond. to (α t , β t , γ t , u t , σ 2 t+1 , k t , l t , q t , y)

r t+1 ∼ G a r + d 2 , b r + 1 2 d ∑ i=1 (k i -q) 2
(iii). Simulate q t+1 cond. to (α t , β t , γ t , u t , σ 2 t+1 , k t , l t , r t+1 , y)

q t+1 ∼ N [σ -2 q + rd] -1 (σ -2 q + r d ∑ i=1 k i ), [σ -2 q + rd] -1
(iv). Simulate l t+1 cond. to (α t , β t , γ t , u t , σ 2 t+1 , k t , q t+1 , r t+1 , y) v). Simulate k t+1 cond. to (α t , β t , γ t , u t , σ 2 t+1 , l t+1 , q t+1 , r t+1 , y) k t+1 ∼ N µ k t+1 , Σ k t+1 (vi). Simulate γ t+1 cond. to (α t , β t , u t , σ 2 t+1 , k t+1 , l t+1 , q t+1 , r t+1 , y)

l t+1 ∼ G a l + d 2 , b l + 1 2 (η t -Kµ H t ) ′ (Σ H ) -1 (η t -Kµ H t ) (
γ t+1 ∼ N µ g t+1 , Σ g t+1
(vii). Simulate β t+1 cond. to (α t , γ t+1 , u t , σ 2 t+1 , k t+1 , l t+1 , q t+1 , r t+1 , y) viii). Simulate α t+1 cond. to (β t+1 , γ t+1 , u t , σ 2 t+1 , k t+1 , l t+1 , q t+1 , r t+1 , y)

β t+1 ∼ N µ b t+1 , Σ b t+1 , B d β + (0, 1) (
α t+1 ∼ N µ a t+1 , Σ a t+1
(ix). Simulate δ t ∼ N (0, Σ MH ), υ t ∼ U [0, 1] and define u t = u t + δ t

• define u t+1 = u t if υ t < π( u t |α t+1 , β t+1 , γ t+1 , σ 2 t+1 , k t+1 , l t+1 , q t+1 , r t+1 , y) π(u t |α t+1 , β t+1 , γ t+1 , σ 2 t+1 , k t+1 , l t+1 , q t+1 , r t+1 , y)

• or u t+1 = u t otherwise where the covariance matrix Σ MH used in this last Metropolis-Hastings step is first estimated over a burn-in phase (the iterations coming from this phase are discarded), and then fixed to its estimated value "asymptotically optimally rescaled" for the final run by a factor 2.38

MCMC algorithm for the estimation of the posterior distribution, using the noninformative prior

In the lines below, we give the different steps of the MCMC algorithm we used to (approximately) simulate (η 1 , . . . , η M ) according to the posterior distribution π(η|y) corresponding to the noninformative prior we presented earlier. Denoting N (µ, Σ, S) the S-truncated Gaussian distribution with mean µ and covariance Σ, the algorithm goes as follows:

Step 1. Initialise η 1 such that π(η 1 |y) = 0

Step 2. For t = 1, . . • or u t+1 = u t otherwise where the covariance matrix Σ MH used in this last Metropolis-Hastings step is first estimated over a burn-in phase (the iterations coming from this phase are discarded), and then fixed to its estimated value "asymptotically optimally rescaled" for the final run by a factor 2.38 d 2 (as recommended for Gaussian proposals in Section 2 of [START_REF] Roberts | Optimal scaling for various Metropolis-Hastings algorithms[END_REF].

The justifications for each full conditional distribution used in the Gibbs sampling steps, including the explicit expressions of µ α t+1 , Σ α t+1 , µ

Figure 1 :Figure 2 :

 12 Figure1: Estimation errors on the four coefficients of the polynomial regression model (9), using MLE (in green), non-informative Bayesian (NI, in blue) and informative Bayesian (I, in red) methods, when ρ = 1 (top row), ρ = 0.5 in the first situation (middle row), ρ = 0.5 in the second situation (bottom row)

Figure 3 :

 3 Figure 3: Estimation errors on the four coefficients of the auto-regressive model (10), using MLE (in green), non-informative Bayesian (NI, in blue) and informative Bayesian (I, in red) methods, for three different sample sizes : N = 100 (top row), N = 30 (middle row), N = 10 (bottom row)

  2. the parameter (a, b) are fixed, we take the values of the estimates calculated on the large sample y H (See Albert (2009); Gelman et al. (2013))

Figure 4 :

 4 Figure 4: Estimation of the hyperparameters θ = (a, b) using hierarchical informative prior (Info.) and non informative prior as function of sample sizes N. The true gamma parameters are (a H , b H ) = (a, b) = (1, 1). The boxplot is constructed from 500 independent replications.

Figure 5 :

 5 Figure 5: Mean square error E n defined in (11) for the three choice of prior : hierarchical prior (Info.), fixed hyperparameter (Fixed) and non-informative prior (ni). The parameters are θ = (a, b): (1, 1) (red), (9/4, 3/2) (green), (4, 2) (blue). The sample sizes are N = 25, 50, 100, 500 and N H = 5000. The gamma parameters of y H are θ H = (a H , b H ) = (1, 1). The boxplot is constructed from 500 independent replications

Figure 6 :

 6 Figure 6: Left : French Electricity load from 13/06/2005 to 29/06/2005 (in grey) and from 05/12/2005 to 11/12/2005 (in black). The load is expressed in MW. Notice the daily patterns of the electricity load are not the same during summer and winter. Right : French Electricity load at 10:00 over 5 years against temperatures. The load seems to increase linearly with the temperature below a certain threshold.

Figure 7 :

 7 Figure7: Quality of the predictions (MAPE in %) averaged for each instant, for populations y and y ′ (left). Posterior densities of the similarity coefficients k j for populations y and y ′ at midday (right).

  proposals in Section 2 of[START_REF] Roberts | Optimal scaling for various Metropolis-Hastings algorithms[END_REF]. The justifications for each full conditional distribution used in the Gibbs sampling steps, including the explicit expressions of µ α t+1 ,

Table 1 :

 1 Estimation y H Estimation y Prediction y Sample size (in days) of the datasets for both experiments.

	Case 1	1099	125		28
	Case 2	1099	144		38
				Case 1		Case 2
			RMSE	MAPE	RMSE	MAPE
	informative prior (I)	770.71	3.08	2041.70	5.71
	non-informative prior (NI)	24440.52 100.26	7317.56	22.03
	non-parametric (NP)	1461.52	5.82	25091.68	77.41
	double exp. smoothing (DES) 1800.12	7.41	22572.58	71.36
	linear time series (ARIMA)	1702.39	6.89	13839.85	43.73

Table 2 :

 2 Overall quality (RMSE in MW, MAPE in %) of the predictions for both experiments.

Table 3 :

 3 Overall quality (RMSE in MW, and MAPE in %) of the estimation (left) and prediction (right) for the non-informative (non-info.) and informative (info.) priors, depending on the number of months used for the estimation (from 12 months to 6 months).

	12 months
	10 months
	8 months
	6 months

  . , M -1, repeat (i). Simulate σ 2 t+1 cond. to (α t , β t , γ t , u t , y) i.e. (ii). Simulate γ t+1 cond. to (α t , β t , u t , σ 2 t+1 , y) i.e. Simulate b t+1 cond. to (α t , γ t+1 , u t , σ 2 t+1 , y) i.e. (iv). Simulate a t+1 cond. to (β t+1 , γ t+1 , u t , σ 2 t+1 , y) i.e. Simulate δ t ∼ N (0, Σ MH ), simulate υ t ∼ U [0, 1] and define u t = u t + δ t • define u t+1 = u t if υ t < π( u t |α t+1 , β t+1 , γ t+1 , σ 2 t+1 , y) π(u t |α t+1 , β t+1 , γ t+1 , σ 2 t+1 , y)

	σ 2 t+1 ∼ IG	N 2	,	1 2	y -f (η) 2 2
	γ t+1 ∼ N µ γ t+1 , Σ	γ t+1
	(iii). β t+1 ∼ N µ	β t+1 , Σ	β t+1 , B + (0, 1) d β
	α t+1 ∼ N µ α t+1 , Σ α t+1
	(v).

Laboratoire de Mathématiques Jean Leray,

Rue de la Houssinière -BP 92208, 44322 Nantes Cedex

3, France 2 Electricité de France R&D, 1 Avenue du Général de Gaulle, 92141 Clamart Cedex, France

Acknowledgments

The authors would like to thank Adélaïde Priou for collecting a part of the data as well as the corresponding results, and Virginie Dordonnat for the insightful discussions.

Appendix

Notations

Using the notation M i• for the i-th row of a matrix M, the non-linear model described in ( 12) can be re-written in the following condensed way, as a particular case of the regression model ( 5): for t = 1, . . . , N,

The matrices

and T of size N × 1 are known exogenous variables while the parameters of the model to be estimated are

where

Technical proofs

Assume furthermore that N > d α + 1 and that (y 1 , . . . , y N ) are observations coming from the model (13) and the posterior measure is then a well-defined (proper) probability distribution.

Proof. First notice that π(η, k, l, q, r|y) dσ 2 is proportional to

for almost every y and that the function θ → yf (θ) -N 2 is bounded, for almost every y. The posterior integrability is hence trivial as long as π(η|k, l)π(k|q, r)π(l)π(q)π(r) itself is a proper distribution which is the case here.


Rodolfo Sergio

González Castro
email: rodolfogonzalezcastro@gmail.com

THEORY OF GRAVITY "ENERGY-WAVE": THE ORIGIN

Keywords: hal.archives-ouvertes.fr:hal-00947254 PACS: 04.50.Kd physics, particles, gravity, gravitation constant, Newton

In this paper I derive an equation relating the gravitational acceleration with the gravitational wavelength corresponding to the "gravitational energy" density at a point in space and with the speed of light, without using the gravitational constant of Newton (G), derive an equation of the Energy-Momentum of Einstein suppressing this constant, and further I set the foundations for new theory of gravity "Energy-Wave".

INTRODUCTION

The mathematical physicist Sir Isaac Newton in 1687 published his book "Philosophiae Naturalis Principia Mathematica" where he presented the law of universal gravitation empirically derived to describe and calculate quantitatively the mutual attraction of each particle and massive objects in the universe. In that document, Newton concluded that the attraction together two bodies is proportional to product of their masses and inversely proportional to the square of the distance that separates them.

However, these must be adjusted proportionalities by introducing a constant called Universal Gravitation (G) with an approximate value of 6.674 × 10 -11 Nm 2 kg -2 units in the International System. Without the introduction of this constant, the equation, lose their rationality and is impractical to calculate the force of gravity without it.

For Newton, gravity is considered a force exerted on a remote Instant. Moreover, when the gravity force is exerted by two or more bodies extremely mass, Newton's law has serious limitations and then must resort to the Theory of Relativity General stated by Albert Einstein in 1915, who says that gravity is not a force exerted to distance but a contraction of Space-Time produced by the presence of Energy-Matter (1).

However in the final formulation of the equations of the universe, to make it compatible with the law of conservation of energy and principles of general covariance, Einstein included geometric concepts such as the Ricci Tensor and scalar, but mainly the Energy-momentum tensor, but fails to integrate into said Energy-Momentum Tensor the constant Universal Gravitation Newton (G), as this is finally out of tensor in the second member of the equation.

While Einstein equation establishes the relationship between gravity, energy and geometry distortions space-time, it does not define the origin of the relationship.

In this regard, in 1995 Jacobson achieved considerable progress in linking the laws of thermodynamics to the Einstein equation and the equation of state correlates entropy with the area of energy flow [START_REF] Jacobson | Thermodynamics of Spacetime: The Einstein Equation of State[END_REF]. Erik Verlinde published on January 6, 2010, his work "On the Origin of Gravity and the Laws of Newton" [START_REF] Verlinde | On the Origin of Gravity and the Laws of Newton[END_REF], which proposed that gravity is a reality entropic force emerging space. In its formulation, includes reduced Planck's constant, N as a Screen of information bits of space, adds a new constant called G, which ultimately found to be equivalent to the Universal Gravitational Constant. On that basis, Verlinde forecast to gravity as a fundamental force.

In March 2010, Jae-Weon Lee, Hyeong-Chan Kim and Lee Jungjai published a paper in which suggest that the Einstein equation can be derived from the Landauer 's principle on the Elimination of information causal horizons, and conclude that gravity has its origin in quantum information [START_REF] Lee | Gravity from Quantum Information[END_REF]. From then such work is also supported by Jacobson linking between thermodynamics and the equation Einstein, as well as on the work o f Verlinde entropic force.

Thus, today we already have a strong linkage between energy, heat, temperature, laws of thermodynamics, general theory of relativity, perturbation of the geometry of space-time, entropy and quantum information, but somehow linking gravity and electrostatic force has failed in all these works that eventually drift to the Einstein equation and repeatedly use the constant of universal gravitation Newton (G).

Derive a gravity equation that eliminates the universal gravitational constant and link the gravitational acceleration with electrostatic acceleration requires first that the equation derived from the principal component of gravity, energy (E), and secondly considers the relation of this energy with Planck area and the emergence of two tension vectors with corresponding emission wavelength that generate electrostatic and gravitational waves as seen in the present work.

Newton's equation establishes the relationship between mass and gravity, while the equation Einstein relates the Energy-Momentum Tensor with the modification or distortion of space:

For particles at absolute rest, Einstein's equation has only one active component ) of the Energy-Momentum Tensor , dimensionally same as is defined by the equation:

(1.1)
Where is the Lorentz factor, c the speed o f light and ƿ the energy density, so if you divide between we simply obtain the energy density, that is, in real terms the Einstein equation defines the energy density is the curving space.

The problem of the Einstein equation is that the energy density in order to be equivalent to the curvature of space, requires the tensor is multiplied by the constant of universa l gravitation Newton (G) and their corresponding dimensions.

To resolve this problem in the Einstein equation, this paper begins from the energy (E) component and distribution or displacement in the area occupied promptly in space under the principle that energy is always a wave packet with a pilot waveform that is proportional to the acceleration of gravity occurs in the space. This will allow us to the end of this study derive an Einstein equation with Energy-Momentum Tensor that eliminate G while promptly defines the relationship between wavelength and space contraction gravitational acceleration.

GRAVITATIONAL ENERGY DERIVED FROM PLANCK AREA

To derive the gravity equation, I need to hypothesize a priori that the particles at absolute rest, energy it's distributed in a specific area of space and that this area has two main components from which emerge two tension vectors on the space with different action and different time.

The first tension vector is considered by Einstein and in principle only generates the electrostatic force and the related electromagnetic phenomena, as the second tension vector becomes the cause of the gravitational interaction. However, the immediate question is where did emerge the second tension vector? To explain the above and get the answer, consider a priori a hypothetical particle whose Plank energy distribution in space is the area of a perfect square:

= = = = 1.03128928E-68 m 2
Both the axis (X) and the axis (Y) corresponding to the Planck length (lp), if we multiply it by 2 then we get the Planck Pilot wavelength ( ). Then immediately notice that the Planck energy density has not one (1) but two (2) waves associated interaction.

Both sides of the square shafts or exert pressure on space, or in terms of Einstein, produce a Tension Vector. In the case of Axis (X) the Vector of Tension is exactly the same as Einstein, but in the case of axis (Y), a hypothetical still correspond Second Tension Vector.

However, in the case of the hypothetical particle of Planck, because either axes or sides of the same length, and in the case of (Y) the Vector Pressure is orthogonal and contractive (from infinity to the X axis), then the two Tension Vectors cancel each other. In terms of the de Broglie wavelength to be exactly equal in length and width but opposite to each other, they cancel each other, so that an external observer cannot feel the presence of the particle of Planck on electrostatic terms or in gravity terms. It is then a particle "Null" or space "Empty".

However, when the axis (X) corresponding to the first Tension Vector (Wavelength De Broglie) extends beyond the Planck length (lp), then proportionally the axis (Y) which corresponds to the Second Tension Vector is shortened (gravitational wave length). But also at the time, in which the wavelengths (X) and (Y) are different, they cease to cancel each other and so the electrostatic and gravitational effects are "visible" to outside observers.

In other words, and under these arguments Planck Area or is always invariant to perform the respective calculations considered here.

For example, let's see here graphically proportional representation of an electron according to Planck area:

Y = = = 1.0155241E-34 m X = = = 1.0155241E-34 m = = 1.03128928E -68 m 2
Note: The graph of the rectangle is not exactly proportional to the scalar quantities for obvious reasons.

As can be observed in the electron because (X) is extended too much relative to the Planck length, then (Y) is shortened proportionately then being a wavelength excessively small but sufficient to generate the gravitational effects of the electron.

In the case of the electron, the De Broglie wavelength is about 2.4263170 E-12 meters, while the gravitational wave length is approximately 4.2504295 E-57 meters.

In the case of the earth, being an object (particle) with a large mass, the proportion is reversed; the De Broglie wavelength is about 3.7009 E-67 meters, while the gravitational wave length is about 0.027866232 meters, according to calculations using the equations, described later.

In other words, the particles and the mass objects have a spatial configuration with respect to the wave packet associated like strings or threads vibrating whose length in the X axis is much longer or shorter with respect to the Planck length but a width or thickness (Y axis) is proportional to the axis X excessively longer or shorter, but enough to be the sources of the gravitational effects..

To understand this, let's draw a conceptually wave example:

The physical interpretation of the Tension Vector I, not a major problem, because it is the product of the tangential acceleration of energy distributed in the wave curve along the wavelength (X axis) with the corresponding pressure on space.

While the physical interpretation of tension Vector II corresponds to the component of normal acceleration (gravitational acceleration) occurring Wave inward, towards the center of mass-energy. In this regard there are three possible interpretations:

I) The energy density distributed in the Wave, has a "thickness or wide" quantum limit set by the Planck Area (yarn thickness wave). In this case, the quantum energy density pressure inside the particle itself. II) A second possible interpretation is that the wave to advance, the vector has a torque space that causes progress turning on itself elliptical. Therefore the rotation or tension energy density causes a second acceleration towards the center of the elliptical rotation, then this acceleration corresponding to the gravitational acceleration. III) In a third interpretation, a second tension vector actually corresponds to the wave amplitude (A) that is not symmetrical with respect to the wavelength, then said amplitude wave then generates a second acceleration whose second component (normal acceleration) is orthogonal to with respect to X (wavelength) axis. Whatever the correct physical interpretation of Tension Vector II, in the end, this also generates a second emission wavelength with a corresponding length.

Moreover to calculate both wavelengths: De Broglie wave length (Tension Vector I), and wave length Gravitational (Tension Vector II), then we have to set the "Energy" which is the source of both.

If we make this definition from the Planck hypothetical particle that transforms into another particle conserving Planck area, then we need to be the axis (X) or Tension Vector I, which corresponds to the energy at absolute rest and which also corresponds to the extension or horizontal length, then we have to multiply the Energy Plank, quantified for this energy density at absolute rest. While the axis (Y) or Tension Vector II corresponds to the "gravitational energy", which in turn corresponds to the width or orthogonal extension, then we also divided quantified according to the following equation:

= (2.1)
1.9560852E+09 Joules

3.82627E+18 Joules 2

Where Ep is Planck's energy and N is a quantum number that corresponds to the ratio of the rest energy of any particle or object mass and the Planck energy given by the equation

(2.2)
Substituting N into Equation (2.1) then the energy distribution in an area Planck is given by:

(2.3)
In this case (E) corresponds to the rest energy of any known particle or massive object, while is the "Gravitational Energy".

However, here I must clarify that the "gravitational energy", which is excessive and inexplicable in the case of the electron and the proton, does not correspond to an existing energy within the particle, because if so would violate the very principles of Plank, but which is the sum of quanta energies of all "Plank particles" which are displaced towards the axis (X) by the contraction of (Y).

By then calculating said sum of displaced quanta energies of Planck (Y) to the axis (X) occupies the following equation:

(

Should clarify here, first, that the boundaries of physics established by the laws of Planck, not defined in linear terms of length or Planck energy, ( , Ep) because then we would have a universe whose total contained energy not would be greater than Ep.

In that sense, the limits of Planck must actually correspond to the Planck Area (Ap) and about the limits of energy; the multiplication of the energy corresponding to the Tension Vector I (E) for the energy corresponding to the Tension Vector II , and whose both product becomes Planck energy squared as defined by equation 2.1 and 2.3.

Moreover, to calculate the "gravitational energy" of an electron, a proton, a neutron or other particle those are located in the "absolute vacuum" we get a magnitude of such "gravitational energy" which apparently violates the limit of Planck (Ep) but not the limit in equation 2.1 or 2.3.

I should also clarify that Tension Vector I corresponding to the energy at absolute rest (E) has a direction from the particle to infinity, so that almost all the energy is contained within the particle itself.

Meanwhile, the Tension Vector II corresponding to the "gravitational energy" has an orthogonal inverse direction from infinity to the particle, so that the majority of the "gravitational energy" is outside of the particle.

Thus "gravitational energy" is the sum of quanta energies of Planck that is moved from infinity to the Tension Vector I (axis X) by the contraction of Tension Vector II (Y) at a rate defined by equation 2.4.

Both the rest energy (E) as the "gravitational energy" (Eg) is "moving" or "distributed" along the corresponding axis (X, Y), and thus generate their own wave packet with pilot wave associated to the front.

The two Pilot waves have a corresponding wavelength according to the principles of Broglie known equation given by: That in the case of "gravitational energy" or Tension vector II, its wavelength can be obtained directly from the following equation:

(2.5)

Where (E) corresponds to the energy at rest.

The product of two wave lengths, corresponding to the wavelength squared Planck. As Planck Area, the product of both wavelengths is invariant in all of the particles:

(2.6)
Where h is Planck's constant and c is the speed of light and where the wavelength of Broglie (Tension Vector I) is and is the gravitational wave length (Tension Vector II).

GRAVITATIONAL COUPLING

Because the two wavelengths are extremely different in most massive particles or objects, as these originate practically in a "simultaneous", the gravitational wave length to be coupled to the length of wave energy , which generates a coupling factor for the gravitational force given by:

(3.1)
Where is the factor of gravitational value coupling.

For the electrostatic force, value coupling is considered constant (Fine Structure) about 7.297352568E-03.

In the case of gravity, said coupling depends on the ratio of both wavelengths as shown in the equation 3.1.

Defining the gravitational engagement is critical to the development of a theory of gravity on the equivalence to approach the electrostatic force and to calculate the force of gravity, but such coupling is not necessary to calculate the gravitational acceleration on the arguments presented here and as it will show at the end of this paper.

THE EQUATION OF GRAVITY FORCE WITHOUT CONSTANT GRAVITATION

In short, the "gravitational energy" ( ) "moves" or "distributed" from infinity to the center of the distribution of particle energy, generating a curvature or contraction of space defined by the length "gravitational wave" , in turn generating component of normal acceleration towards the center of the particle which we interpret as gravitational acceleration.

Because of this, we then multiply the "gravitational energy" corresponding gravitational wave length: But because of the wide differential between rest energy and "energy gravitational" or the wide differential between De Broglie wavelength and gravitational wave length, then we have to include in the equation the gravitational coupling factor obtained from equation 3.1:

And then for the gravitational attraction between two identical particles or massive objects or with the same amount of mass only remains for us to divide the radius separating the two particles:

= (4.1)
Where is the "gravitational energy" (equation 2.4), is the length of gravitational wave (equation 2.5), is the gravitational coupling (equation 3.1) and (d) is the radius between the two particles interacting or massive objects.

And for the case of two different particles or objects mass:

= (4.2)
Then comes to be an equation for calculating the gravitational interaction without the use of the gravitational constant of Newton and use the length of de Broglie wave, for particle or massive objects at rest.

But due to the existence of a constant in this equation, it is possible to reduce it further.

1.986451698E-25 joules per meter (4.3) Due to this constant, we can use in the equation 4.2, only the values of "gravitational energy" and gravitational wave length of the particle or object with higher energy (or mass) and then the equation is abbreviated :

4) )
Keeping 2 in the equation is to finally remind us that the ratio between two particles or massive objects is not straight but curved.

Moreover, also in the case of the gravitational interaction between a massive object and a smaller, we can eliminate the two gravitational couplings of equation 4.4 by exchanging the component of the gravitational energy ( for the energy in rest (E) of the massive object while retaining the gravitational wave length :

4.6)
Where (F) is the force of gravitational interaction, (E) the rest energy of the massive object, and the length of pilot gravitational wave associated with the gravitational energy gained by the equation 2.5.

DERIVATION OF THE EQUATION OF GRAVITATIONAL ACCELERATION

Then we can use Equation 4.6 to derive a new equation for calculating the gravitational acceleration of a massive body like the earth.

According to Newton, force (F) is equal to mass (m) by acceleration (a):

=

We make the Newton equation equal to equation 4.6:

(5.1)

Convert (E) to the Einstein equation:

(5.2)
Spent the terms:

(5. 3)
Then we exchanged the term corresponding to the mass to remove it from the equation:

(5.4) Eliminates m:

(5.5)

We returned the radius And Eureka:

= (5.6)
That is, the gravitational acceleration is equal to the length gravitational wave of the mass object (eg earth) by the speed of light squared divided radius squared.

As we can see in equation 5.6, noticeably disappear two components: force (F) and mass (m). The immediate interpretation of this is that the energy density and the corresponding wavelength do not exert a force on objects that attracts, but this force is exerted actually about the same surrounding space which is accelerated in an orthogonal contraction towards the center of Energy density in proportion to the length of dominant gravitational wave. Substituting the values for the gravitational acceleration of the earth:

Ensure equation in the

= =
Which is the same result obtained with Newton's equation.

DERIVATION UNIVERSAL GRAVITATIONAL CONSTANT (G)

Conventionally the gravitational acceleration is calculated with:

=

We then said equation equivalent to equation 5.6:

= (6.1)
Eliminated radius from the two equations :

= (6.2)
And finally solve for G:

(6.3)
We can convert the units of Planck equation and result is the same:

(6.4)

Where (

) is the Planck wave length defined in 2 and is equivalent to the mass Planck of 2.1764383 E-08 kg.

If we calculate with different measures o f mass or energy using equation 6.3, then we obtain an integer value of G, ie a purely quantum value:

DERIVATION OF EQUATION CALCULATING ELECTROSTATIC INTERACTION

Over the original arguments, calculating the electrostatic interaction between two particles is even simpler because the constant of electrostatic coupling or fine structure. In this case leads us to a broad general constant that shall term Universal Electrostatic Constant:

1.4495849660E-27 joules per meter (7.1)

No matter how much rest energy of the particle, the electrostatic force is always the same and only varies depending on the distance. In this case, this energy corresponds to the Tension vector I (X axis of Planck area)

That is, the calculation of the electrostatic interaction between two identical or different particles is given by the general equation:

= = (7.2)
Where E is the energy at rest; the wavelength of de Broglie; coupling constant or Fine structure; d is the radius between the two particles. This equation is equivalent to: = Thus Coulomb's law for electrostatic interaction between only two particles (electron-electron; proton-electron; protonproton). we have the following equivalence: = =

DERIVATION OF EQUATION OF ELECTROSTATIC ACCELERATION

Following the same reas oning as in 5, then we have:

(8. That is, the electrostatic acceleration equals the de Broglie wavelength of the particle by the speed of light between the squared distances between the interacting particles.

Since equations can be observed 8.6 and 5.6, are virtually identical, with the slight difference that in the case of the electrostatic acceleration have to stop for the time the fine structure constant.

Despite the existence of the fine structure constant in equation 8.6, we can also say that the electrostatic force is not exerted on the mass of the second particle interaction but on its surrounding space which is accelerated so repulsive or attractive depending of the signs of the charge of the particles in interaction.

GRAVITATIONAL ACCELERATION EQUIVALENT TO ELECTROSTATIC ACCELERATION.

The reason why the equations 5.6 and 8.6 are not fully equivalent has a very simple reason; 5.6 is a classic equation, while the second (8.6) is an equation within the framework of the Special Theory of Relativity.

In the equation 8.6, (Fine Structure) actually corresponds to the Lorentz factor of relativistic velocities that cause the contraction of the wavelength, according to the Theory of Special Relativity.

According to Sommerfeld calculations, calculating the velocity of the electron in the first orbit is given by the equation:

= =

Where ( ) is the electrostatic acceleration; ( ) is the gravitational acceleration and ( ) is the Lorentz contraction factor, we find that the only noticeable difference between the two equations is the wavelength.

To derive the equality between the two types of acceleration, we only have to use N, the quantum number obtained in equation 2.2: = (9.7)

DERIVATION OF EQUATION OF EINSTEIN TENSOR

Then only remains for us is to derive an equation of Einstein Tensor that integrates the arguments presented here and that otherwise removed from the equation the gravitational constant as posed in the introduction to this work.

We use the active component of the Einstein tensor in the case of particles at rest or in nonrelativistic velocities:

We took that component of the energymomentum tensor: But now it has been the equation in terms of the wavelength of the particle or massive object that contains all the information practically the same:

(10.9)

And if we want to integrate the space-time coordinates:

(10.10)

Where (A) is the area occupied by the particle, (c) the speed of light and (t) is the time.

Finally we have the modified equation of Energy-momentum Tensor of Einstein and gravitational constant of Newton deleted from it:

CONCLUSION

Having derived clean and natural equation gravity and shape of the same equation gravitational acceleration without using the gravitational constant of Newton, with a different way to Newton and Einstein, using as components only terms of energy and wavelength; or speed of light and wavelength, the way for a new interpretation of gravity opens, even for a new theory of gravity "energy-wave":

1. That to derive the equations of gravity in terms electrostatic necessarily need to include the concept of "Gravitational Energy" and the concept of Second Tension Vector derived Planck Area and the limit defined by the Planck Energy Square.

2. That it is possible to derive equations for the electrostatic acceleration and gravitational acceleration on the basis of the same principles.

3. That in the case of the force equation Gravity, because the length of wave pilot of resting energy and "gravitational energy" wave pilot length, are vastly different and also originate nearly so "simultaneous"; "Gravitational" wave length should fit variable length electrostatic wave form, which had hitherto hindered its calculation electrostatically.

4. However, to derive the equation of gravitational acceleration (5.6), gravitational coupling factor disappears, generating a surprising equivalence of gravitational wave length and gravitational acceleration. 5. Energy density and its corresponding wavelength does not exert a gravitational force on objects that attracts, but this force is exerted on it actually surrounding space is so rapid orthogonal contraction towards the center of density of energy in proportion overriding the gravitational wave length.

6. The electrostatic force is not exerted on the mass of the second particle in interaction but surrounding space there on which is accelerated in a repulsive or contraction depending on the signs of the charge of the particles interact.

7. That the final equation in 8.6 here developed for the estimation of the gravitational acceleration, is actually an electrostatic acceleration equation where the fine structure constant is the Lorentz contraction factor for the wavelength in the context of the special theory Relativity of Einstein.

  Squared radius of the Earth: 4.05896E+13 m

= (9. 1 )= 3 )(9. 4 )

 134 Where (h) is Planck's constant. Substituting the values: = We now apply the Lorentz factor in summary, because here we have obtained a speed within the framework of the theory of special relativity:That becomes the Lorentz contraction factor for the wavelength of an electron the first Bohr orbit, and is constant value for successive orbits or energy levels of electrons in the atom. So in this case:In basing the above arguments then modify equation 8.6 for electrostatic acceleration:=Then the equation of the gravitational acceleration on the 5.6, also transform to the framework of the Special Theory of Relativity including the Lorentz contraction factor for gravitational wave length: ), is the traditional symbol of the gravitational acceleration.And finally we have the equivalence:

  wave length obtained in equation 2.5 and ( ) is the Lorentz contraction factor for relativistic particles, and (V) is the volume occupied by the mass-energy in space.