Wilfried Kirschenmann
email: wilfried.kirschenmann@edf.fr

Laurent Plagne
email: laurent.plagne@edf.fr

Stéphane Vialle
email: stephane.vialle@supelec.fr

Multi-Target Vectorization With MTPS C++ Generic Library

Keywords: GPU, SSE, Vectorization, c++ Template Metaprogramming, Performances

This article introduces MTPS, a C++ template library dedicated at vectorizing algorithms for different target architectures. Algorithms written with MTPS benefit from optimized memory access patterns and show performances close to hardware limits, both on multicore CPU and on GPU.

Introduction

In many scientific applications, computation time is a strong constraint. Optimizing these applications for the rapidly changing computer hardware is a very expensive and time consuming task. Emerging hybrid architectures tend to make this process even more complex.

The classical way to ease this optimization process is to build applications on top of High Performance Computing (HPC) libraries that are available on a large variety of hardware architectures. Such scientific applications, whose computing time is mostly consumed within such HPC library subroutines, then automatically exhibit optimal performances for various hardware architectures.

However, most classical HPC libraries implement fixed APIs (e.g., BLAS) and may be too rigid to match the needs of all client applications. In particular, classical APIs are limited to manipulate rather simple data structures like dense linear algebra matrices. As a more complex issue, general sparse matrices cannot be represented with a unified data structure and various formats are proposed by more specialized libraries. In the extreme case, structured sparse matrices cannot be efficiently captured by any of the classical library data structures. Relying on such complex matrices, several neutron transport codes developed at EDF R&D require another kind of library to be used.

Following the model of the C++ Standard Template Library (STL), template based generic libraries such as Blitz++ [START_REF] Veldhuizen | Arrays in blitz++[END_REF] provide more flexible APIs and extend the scope of library-based design for scientific applications.

Such generic libraries allow to define Domain Specific

Diagonal Tridiagonal

Examples of composition Embedded Languages (DSELs) [START_REF] Czarnecki | Dsl implementation in metaocaml, template haskell, and c++[END_REF]. Legolas++, a basis for several HPC codes at EDF, is a C++ DSEL dedicated to structured sparse linear algebra. In order to meet EDF's industrial quality standards, a multi-target version of Legolas++, currently under development, will provide a unified interface for the different target architectures available at EDF, including clusters of heterogeneous nodes (i.e., with both multi-core CPUs and GPUs). This article presents MTPS (Multi-Target Parallel Skeletons), a C++ generic library dedicated to multi-target vectorization that is used to write the multi-target version of Legolas++. Only developments concerning a single heterogeneous node are presented here.

The next section presents the principles of Legolas++ and Section 3 introduces MTPS. Its optimization strategies and the achieved performances are discussed in Section 4. Finally, conclusions are drawn in Section 5.

Towards a Multi-Target Linear Algebra Library

Legolas++ is a C++ DSEL developed at EDF R&D to build structured sparse linear algebra solvers. Legolas++ provides building bricks to describe structured sparse matrix patterns and the associated vectors and algorithms.

Legolas++ is based on the observation that most structured sparse matrix patterns can be described as the composition of simpler patterns. Figure 1 shows how to compose two simple patterns to create new patterns. Figure 2 shows an example of matrix pattern issued from a neutron transport code written with Legolas++ (GLASS in [START_REF] Plagne | Generic programming for deterministic neutron transport codes[END_REF]).

The explicit GPU parallelization of one of our neutron transport code resulted in speed-ups around 30 over the sequential Legolas++ CPU implementation [START_REF] Kirschenmann | Massively parallel solving of 3D simplified P N equations on graphic processing units[END_REF]. To generalize this gain of performances to other Legolas++ based applications, a parallel and multi-target version of Legolas++ is being developed. As the parallel CPU and GPU versions exhibit strong similarities, Legolas++ developments Some libraries, like TrilinosNode [START_REF] Baker | A light-weight api for portable multicore programming[END_REF], Quaff [START_REF] Falcou | Quaff: efficient c++ design for parallel skeletons[END_REF] or Intel TBB [START_REF] Reinders | Intel threading building blocks[END_REF], require their users to explicitely express the parallelism within the application by using parallel skeletons. This expression of available parallelism can be encapsulated into specialized and implicitely parallel STL-like containers and algorithms, as in Thrust 1 and Honei [START_REF] Van Dyk | Honei: A collection of libraries for numerical computations targeting multiple processor architectures[END_REF].

Our goal is to provide implicit parallelism within Lego-las++ containers and algorithms. To ease the writing of its containers and algorithms, Legolas++ relies on MTPS which follows a parallel skeletons based approach. Then MTPS optimizes the code for the different architectures.

As this article presents MTPS, only code for MTPS is shown. For Legolas++ users, MTPS details are hidden in its containers and algorithms.

Collections and Vectorizable Algorithms

This section introduces the notions of collection and vectorizable algorithm on which MTPS relies.

In C++, a Plain Old Data (POD) is a type whose memory representation pattern can be changed without altering its value [START_REF]Programming languages -C++. International Organization for Standardization[END_REF]. POD members can be either integral types or PODs. In the following code snippet, MyPOD is a POD with three float data members:

1 struct MyPOD { float a , b , c ; } ;
Let a collection be a data structure containing different instances of the same POD and f be a pure function (i.e., f has no side effects). An algorithm applying f to all elements of a collection is said to be vectorizable. To parallelize such algorithms, MTPS provides two parallel skeletons optimized for different target architectures: map and fold which correspond to a parallel for loop and to a parallel reduction respectively.

An algorithm is vectorizable in reference to a given collection. For instance, an algorithm operating on a matrix can be either vectorizable or not depending on whether the matrix is considered as a collection of columns or as a collection of rows. Two algorithms vectorizable in reference to the same collection are said to be in the same vectorial context. On the contrary, if two consecutive algorithms are not vectorizable in reference to the same collection, a context switch is required. In a distributed memory system, context switches correspond to communications.

Linear Algebra Hello World of MTPS: saxpy

This section presents the MTPS implementation for the saxpy operation whose implementation in C is:

1 float * X , * Y , a ; 2 for ( int i = 0 ; i<N ; ++i ) Y [ i ]+= a * X [ i ] ;
First, the iteration-dependent data are gathered in a POD XYData whose members correspond to

X[i] and Y[i].
The types of the two members (float) are passed as template arguments to MTPS::POD and their names (x and y) are given in the Fields enum: Second, a collection of XYData elements, xyCol, can be built using MTPS containers. A class type for optimized container is provided as member of the class corresponding to the target architecture. Two levels of parallelism are available on CPUs: thread parallelism and SIMD parallelism. The choice for each level is made by passing two arguments to the CPU template class. Thread can be one of MTPS::Sequential, MTPS::OMP (openMP) or MTPS::TBB (Intel TBB). SIMD can be one of MTPS::Scalar or MTPS::SSE. On CUDA-enabled GPUs, only the SIMD parallelism is used.

1 / / t y p e d e f MTPS : : GPU : : CUDA T a r g e t ; 2 typedef MTPS : : CPU<Thread , SIMD> Target ; Third, the function that is to be applied to all elements of the collection must be written as a functor class AxpyOp which contains the value a internally: As XYData elements may not be stored identically on different target architecture, AxpyOp::operator() does not take an XYData as argument. A View is provided instead. XYData members can be accessed with the operator() of the View which takes an int as argument. Elements of the Fields enum can be used either to initialize an int (line 7) or directly (line 8). The declaration of AxpyOp::operator() must be preceeded by the INLINE macro which defines target-dependent keywords (e.g. device for CUDA).

Finally, the functor can be passed to the map and fold parallel skeletons provided by the collection container:

Optimization of performances

For each architecture, the specific optimizations required to enable good performances will be presented. The implementation of an example will then be discussed.

Multi-Target Performance Optimizations

Parallelizing a vectorizable algorithm is straightforward. However, achieving good performances is not: modifications of the collection storage pattern may be required. Indeed, achieving efficient usage of memory bandwidth on a given hardware architecture requires specific access patterns [START_REF] Kirschenmann | Multi-target c++ implementation of parallel skeletons[END_REF]. Figure 4 shows how a matrix of 8 TriDiagonal Symmetric (TDS) blocks of size 4 is stored on three different architectures to optimize the memory access pattern.

Performances achieved thanks to this optimization will be shown in Section 4.3. As this optimization is made in MTPS collection container, MTPS user must define both the size per POD-element of each field (4 for the diagonal field on Figure 4) and the number of POD-element in order to construct a collection. Using these two information, MTPS optimizes the storage for each target architecture.

A context switch imply a data reordering. For instance, switching a collection of matrix rows to a collection of matrix columns modifies the effective storage (i.e. the matrix is transposed). This part of MTPS is under development.

Implementation of a Linear System Resolution

The example presented in this section corresponds to a basic operation that represents the major part of the execution time of a neutron transport code [START_REF] Kirschenmann | Massively parallel solving of 3D simplified P N equations on graphic processing units[END_REF]. Let A be a blockdiagonal matrix with TDS blocks. smaller problems. The AX = B linear system can be seen as a collection of smaller block systems ax = b that can be solved independently. To solve one ax = b system, the matrix a is factorized in-place with a LDL T factorization and a forward and backward substitution is then applied on x. Only the code for the factorization is shown here.

Let us introduce TData which represents a TDS block. TData elements are stored in two vectors corresponding to the diagonal and the lower diagonal: The Shape type of line 5 contains the effective sizes of the two fields. All elements of a collection have the same shape. As both the number of TData elements and their shape is known, tCol storage pattern can be optimaly stored according to the target architecture (see Figure 4):

1 TData : : Shape s=TData : : createShape ( size ) ; 2 Target : : collection<TData> tCol ( N , s ) ;

The following code snippet corresponds to TLDLtOp which factorize the matrix a in-place using a LDL T decomposition:

1 struct TLDLtOp { 2 template <template <class> class View> The elements of a field can be accessed by passing their index as the second argument of the view operator(). If this index is not provided, its default value is 0. Line 7 shows how the type of a field elements can be retrieved.

int i = 1 ; i < size ; i ++) { 10 l=a ( low , i-1) / a ( diag , i-1) ; 11 a ( low , i-1)=l ; 12 a ( diag , i )-=a ( diag , i-

Performances

Table 1 shows the performances obtained to solve the AX = B system with A having 10 5 blocks of size 100. Speed-ups are given compared to the sequential scalar CPU version. CPU tests are run on a machine with two 2 GHz Intel E5504 quad-core processors. GPU tests are run on a Nvidia Quadro FX5800 card. Both architectures were launched at the end of 2008. Computation performances are given in GFlops. Data throughput is given in GB/s and takes into account the data transfers to and from the memory: on CPU, an element remaining in cache between two loads is considered to have been loaded only once. The achieved performances are compared to the measured peak performances. Peak computational power is measured with large BLAS matrix-matrix multiplications (sgemm): 55 GFlops on CPU and 348 GFlops on GPU. Peak memory throughputs are measured with the stream benchmark [START_REF] Mccalpin | Memory bandwidth and machine balance in current high performance computers[END_REF] on CPU (24 GB/s) and with the CUDA SDK bandwidth benchmark on GPU (74 GB/s).

On both CPU and GPU, parallel performances are limited by the memory bandwidth and our approach reaches more than 95% of the peak memory bandwidth utilization.

Figure 1 :

 1 Figure 1: Matrix pattern composition within Legolas++.

Figure 2 :

 2 Figure 2: A matrix pattern built with Legolas++.

3

  Introduction to MTPS 3.1 Related Work Many libraries parallelize for different architectures from a single source code. A complete bibliography is beyond the scope of this abstract; only some examples based on C++ meta-programming techniques are discussed.

1Figure 3 :

 3 Figure 3: Our hourglass software architecture to achieve a multi-target Legolas++: a minimal MTPS library adapts the code for different hardware architectures.

1 struct XYData 2 :

 12 public MTPS : : POD<float , float>{ 3 enum Fields {x , y } ; 4 } ;

3 4

 3 Target : : collection<XYData> xyCol ( N ) ;

Figure 4 :

 4 Figure 4: The storage of the diagonal is adapted by MTPSfor the target architecture.

1 struct TData 2 : 6 7 8 Shape out ; 9 out

 12689 public MTPS : : POD<float , float>{ 3 enum Fields {diag , low } ; 4 typedef MTPS : : POD<float , float> Base ; 5 typedef typename Base : : Shape Shape ; static Shape createShape ( int size ) { [ diag ] = size ; 10 out [ low ] = size -1; 11 } 12 } ;

3 INLINE 5 typedef

 35 void operator ( ) ( View<TData> a ) 4 const { View<TData> TV ; 6 int low = TV : : low , diag = TV : : diag ; 7 typename TV : : template Type<low > : : Type l ; 8 int size = a . shape ( ) [ diag ] ; 9 for (

Table 1 :

 1 Performances of MTPS for the TDS example. Computation are carried out in single precision floating point. The smallest time over 1000 executions is given.

	1) * l * l ;

Conclusions and perspectives

We have presented MTPS, a C++ generic library simplifying the parallelization and the optimization of vectorizable algorithms for different architectures. In particular, one can program an algorithm once, compile it for execution on the SSE units of multicore CPUs or on CUDA-enabled GPUs and obtain performances close to hardware limits: 95% of peak performances were achieved.

For further developments of MTPS, two main directions are considered. On the one hand, having an efficient implementation of context switches between different vectorizable algorithms is essential to write complex applications. On the other hand, targeting other architectures, including distributed memory architectures, will add more opportunities to increase the performances. These two directions of research will be investigated in a near future.


Pierre Rouchon
email: pierre.rouchon@mines-paristech.fr

May

On stability of discrete-time quantum filters

Fidelity is known to increase through a Kraus map: the fidelity between two density matrices is less than the fidelity between their images via a Kraus map. We prove here that, in average, the square of the fidelity is also increasing for a quantum filter: the square of the fidelity between the density matrix of the underlying Markov chain and the density matrix of its associated quantum filter is a super-martingale. Thus discrete-time quantum filters are stable processes and tend to forget their initial conditions.

Kraus maps and quantum Markov chains

Take the Hilbert space S = C n of dimension n > 0 and consider a quantum channel described by the Kraus map (see [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF], chapter 4)

K(ρ) = m µ=1 M µ ρM † µ (1)
where • ρ is the density matrix describing the input quantum state, K(ρ) being then the output quantum state; ρ ∈ C n×n is a density matrix, i.e., an Hermitian matrix semipositive definite and of trace one;

• for each µ ∈ {1, . . . , m}, M µ ∈ C n×n /{0}, and µ M † µ M µ = I.
To this quantum channel is associated the following discrete-time Markov chain:

ρ k+1 = M µ k (ρ k ) ( 2 )
where

• ρ k is the quantum state at sampling time t k and k the sampling index (t k < t k+1 ).

• µ k ∈ {1, . . . , m} is a random variable;

µ k = µ with probability p µ (ρ k ) = Tr M µ ρ k M † µ . 1 Tr(MµρM † µ) M µ ρM † µ = 1 pµ(ρ) M µ ρM † µ .
Kraus maps are contractions for the trace distance, i.e., nuclear distance (see [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF], theorem 9.2, page 406): for all density matrices σ, ρ, one has

Tr (|K(σ) -K(ρ)|) ≤ Tr (|σ -ρ|)
where, for any Hermitian matrix A with spectrum {λ l } l∈{1,...,n} , Tr (|A|) = n l=1 |λ l |. The Kraus map tends also to increase fidelity F (see [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF], theorem 9.6, page 414): for all density matrices ρ and σ, one has

Tr K(σ)K(ρ) K(σ) = F (K(σ), K(ρ)) ≥ F (σ, ρ) = Tr √ σρ √ σ (3)
where, for any Hermitian semi-positive matrix A = UΛU † , U unitary matrix and Λ = diag{λ l } l∈{1,...,n} ,

√ A = U √ ΛU † with √ Λ = diag{ √ λ l } l∈{1,.
..,n} . The conditional expectation of ρ k+1 knowing ρ k is given by the Kraus map: , = K(ρ). In section 2, we show during the proof of theorem (1) the following inequality

E (ρ k+1 /ρ k ) = K(ρ k ).
m µ=1 Tr M µ ρM † µ F 2 MµσM † µ Tr(MµσM † µ)
,

MµρM † µ Tr(MµρM † µ) ≥ F 2 (σ, ρ) ( 4 )
for any density matrices ρ and σ. The left-hand side is related to a conditional expectation. Inequality (4), attached to the probabilistic mapping (2), can be seen as the stochastic counter-part of inequality (3) attached to the deterministic mapping [START_REF] Haroche | Exploring the Quantum: Atoms, Cavities and Photons[END_REF]. When for some µ, Tr M µ σM † µ = 0 with Tr M µ ρM † µ > 0, one term in the sum at the left-hand side of (4) is not defined. This is not problematic, since in this case, if we replace

MµσM † µ Tr(MµσM † µ) by MµξM † µ Tr(MµξM † µ)
where ξ is any density matrix such that Tr M µ ξM † µ > 0, this term is then well defined (in a multi-valued way) and inequality (4) remains satisfied for any such ξ.

During the proof of theorem (8), we extend this inequality to any partition of {1, . . . , m} into p ≥ 1 sub-sets P ν :

ν=p ν=1 Tr µ∈Pν M µ ρM † µ F 2 µ∈Pν MµσM † µ Tr( µ∈Pν MµσM † µ) , µ∈Pν MµρM † µ Tr( µ∈Pν MµρM † µ) ≥ F 2 (σ, ρ) (5)
2 The standard case.

Take a realization of the Markov chain associated to the Kraus map K. Assume that we detect, for each k, the jump µ k but that we do not know the initial state ρ 0 . The objective is to propose at sampling k, an estimation ρk of ρ k based on the past detections µ 0 , . . . , µ k-1 . The simplest method consists in starting from an initial estimation ρ0 and at each sampling step to jump according to the detection. This leads to the following estimation scheme known as a quantum filter:

ρk+1 = M µ k (ρ k ) ( 6 ) with p µ (ρ k ) = Tr (M µ ρ k M µ ) as probability of µ k = µ. Notice that when Tr (M µ k ρk M µ k ) = 0, M µ k (ρ k )
is not defined and should be replaced by M µ k (ξ) where ξ is any density matrix such that Tr M µ k ξM µ k > 0 (take, e.g., ξ = 1 n I d ). The theorem here below is a first step to investigate the convergence of ρk towards ρ k as k increases.

Theorem 1. Consider the Markov chain of state (ρ k , ρk ) satisfying (2) and (6). Then

F 2 (ρ k , ρ k ) is a super-martingale: E (F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k )) ≥ F 2 (ρ k , ρ k ).
When ρk or ρ k are pure states, ρk+1 or ρ k+1 remain also a pure states. Then, F 2 (ρ k , ρ k ) = Tr (ρ k ρ k ) and F 2 (ρ k+1 , ρ k+1 ) = Tr (ρ k+1 , ρ k+1 ). In this case, theorem 1 has been proved in [START_REF] Mirrahimi | Feedback generation of quantum Fock states by discrete QND measures[END_REF] using Cauchy-Schwartz inequalities for m = 2. The proof proposed here below deals with the general case when both ρ k and ρk can be mixed states. It relies on arguments similar to those used for the proof of theorem 9.6 in [START_REF] Nielsen | Quantum Computation and Quantum Information[END_REF].

Proof. ρ and ρ are associated to the Hilbert space S = C n : ρ and ρ are operators from S to S. Take a copy Q = C n of S and consider the composite system living on S⊗Q ≡ C n

F (ρ k , ρ k ) = | ψk |ψ k |. We have E F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k ) = m µ=1 p µ (ρ k )F 2 (M µ (ρ k ), M µ (ρ k )).
The matrices M µ (ρ k ) and M µ (ρ k ) are also density matrices. Take the space S ⊗ Q ⊗ E where E is the Hilbert space of the environment appearing in the system-environment model of the Kraus map (1). This model is recalled in appendix A. It introduced an unitary transformation U on S ⊗ E. This unitary transformation can be extended to

S ⊗ Q ⊗ E ≡ S ⊗ E ⊗ Q by setting V = U ⊗ I (I is identity on Q). Then M µ (ρ k ) = Tr Q⊗E P µ V (|ψ k ψ k | ⊗ |e 0 e 0 |) V † P µ . Set |φ k = |ψ k ⊗ |e 0 ∈ S ⊗ Q ⊗ E and |χ k = V |φ k . Using P 2 µ = P µ , we have p µ (ρ k ) = Tr (M µ (ρ k )) = φ k |V † P µ V |φ k = P µ |χ k 2 For each µ, the state |χ kµ = 1 √ pµ(ρ k ) P µ |χ k is a purification of M µ (ρ k ): M µ (ρ k ) = Tr Q⊗E (|χ kµ χ kµ |) . Similarly set | φk = | ψk ⊗ |e 0 and | χk = V | φk . For each µ, | χkµ = 1 √ pµ(ρ k ) P µ | χk is also a purification of M µ (ρ k ). By Uhlmann's theorem, F 2 (M µ (ρ k ), M µ (ρ k )) ≥ | χkµ |χ kµ | 2 .
Thus we have

E F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k ) ≥ m µ=1 p µ (ρ k ) | χkµ |χ kµ | 2 . Since V is unitary, | χk |χ k | 2 = | φk |φ k | 2 = | ψk |ψ k | 2 = F 2 (ρ k , ρ k ). Let us show that m µ=1 p µ (ρ k ) | χkµ |χ kµ | 2 ≥ | χk |χ k | 2 . We have p µ (ρ k ) | χkµ |χ kµ | 2 = | χkµ |P µ χ k | 2 = | χkµ |χ k | 2 ,
| χk |χ k | 2 = | χk | P |χ k | 2 ≤ χk 2 P |χ k 2 = P |χ k 2 = m µ=1 | χkµ |χ k | 2 .
3 The aggregated case.

Let us consider another Markov chain attached to the same Kraus map (1) and associated to a partition of {1, . . . , m} into p ≥ 1 sub-sets P ν (aggregation of several quantum jumps via "partial Kraus maps"):

ρ k+1 = 1 Tr µ∈Pν k Mµρ k M † µ   µ∈Pν k M µ ρ k M † µ   (7)
where ν k = ν with probability Tr

µ∈Pν M µ ρ k M † µ . Consider the associated quantum filter ρk+1 = 1 Tr µ∈Pν k Mµ ρk M † µ   µ∈Pν k M µ ρk M † µ   (8)
where the jump index ν k coincides with the jump index ν k in (7). Then we have the following theorem.

Theorem 2. Consider the Markov chain of state (ρ k , ρk ) satisfying (7) and (8). Then

F 2 (ρ k , ρ k ) is a super-martingale: E (F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k )) ≥ F 2 (ρ k , ρ k ).
Proof. It is similar to the proof of theorem 1. We will just point out here the main changes using the same notations. We start from

E F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k ) = p ν=1 pν (ρ k )F 2 ( Mν (ρ k ), Mν (ρ k )).
where we have set

pν (ρ) = Tr µ∈Pν M µ ρM † µ , Mν (ρ) = 1 pν (ρ) µ∈Pν M µ ρM † µ .
With Pν the orthogonal projector on S⊗Q⊗span{|µ , µ ∈ P

  This result from the trivial identity m µ=1 Tr M µ ρM †

  2 . Then ρ and ρ correspond to partial traces versus Q of projectors | ψ ψ| and |ψ ψ| associated to pure states | ψ and |ψ ∈ S ⊗ Q: ρ = Tr Q | ψ ψ| , ρ = Tr Q (|ψ ψ|) | ψ and |ψ are called purifications of ρ and ρ. They are not unique but one can always choose them such that F (ρ, ρ) = | ψ|ψ | (Uhlmann's theorem). Denote by | ψk and |ψ k such purifications of ρk and ρ k satisfying

  thus it is enough to prove that m µ=1 | χkµ |χ k | 2 ≥ | χk |χ k | 2 . Denote by R ⊂ S ⊗ Q ⊗ E the vector space spanned by the ortho-normal basis (| χkµ ) µ∈{1,...,m} and by P the projector on R. Since | χk = m µ=1 P µ | χk = m µ=1 p µ (ρ k )| χkµ | χk belongs to R and thus | χk |χ k | 2 = | χk | P |χ k | 2 . We conclude by Cauchy-Schwartz inequality

2 For 1 √Similarly | χkν = 1 √F 2 ( 2 . 2 .

 211222 ν } and Mν (ρ) = µ∈Pν M µ ρM † µ , we haveMν (ρ k ) = Tr Q⊗E Pν V (|ψ k ψ k | ⊗ |e 0 e 0 |) V † Pν and pν (ρ k ) = Tr Mν (ρ k ) = φ k |V † Pν V |φ k = Pν |χ k each ν, the state | χkν = pν (ρ k ) Pν |χ k is a purification of Mν (ρ k ): Mν (ρ k ) = Tr Q⊗E (| χkν χkν |) . pν (ρ k ) Pν | χk is also a purification of Mν (ρ k ). By Uhlmann's theorem, Mν (ρ k ), Mν (ρ k )) ≥ | χkν | χkν | Thus we have E F 2 (ρ k+1 , ρ k+1 )/(ρ k , ρ k ) ≥ p ν=1 pν (ρ k ) | χkµ | χkµ | Let us show that p ν=1 pν (ρ k ) | χkν | χkν | 2 ≥ | χk |χ k | 2 = F 2 (ρ k , ρ k ). We have pν (ρ k ) | χkν | χkν | 2 = | χkν | Pν χ k | 2 = | χkν |χ k | 2 ,

Concluding remarks

Theorems 1 and 2 are still valid if the Kraus operators M µ depend on k. In particular, F ( ρk , ρ k ) remains a super-martingale even if the Kraus operators depend on ρk , i.e., in case of feedback. When σ and ρ are pure states (projectors of rank one), D(σ, ρ) = 1 -F 2 (σ, ρ). Consequently inequality (4) yields to

for any pure states σ and ρ (use the fact that [0, x] x → √ 1 -x is decreasing and concave). We conjecture that such inequality hold also true for any mixed states and that

A System-environment model

The quantum channel associated to the Kraus map (1) or the Markov chain (2) admits a system-environment model (see [START_REF] Haroche | Exploring the Quantum: Atoms, Cavities and Photons[END_REF], chapter 4 entitled "The environment is watching"). Take the Hilbert space E = C m associated to the environment and the composite system living on S ⊗ E. Take a pure state |φ k ∈ S and its density matrix ρ k = |φ k φ k |. Assume that before detection µ k at step k, the composite system admits the pure state |φ k ⊗ |e 0 where |e 0 is an environment pure state. Take m states |µ forming an orthogonal base of E. Then exists a unitary transformation U (not unique) of S ⊗E such that, for all |φ ∈ S,


Jean-Bernard Chatelain
email: jean-bernard.chatelain@univ-paris1.fr

Kirsten Ralf

Inference on Time-Invariant Variables using Panel Data: A Pre-Test Estimator with an Application to the Returns to Schooling

Keywords: JEL classification numbers: C01, C22, C23 Time-Invariant Variables, Panel data, Time-Series Cross-Sections, Pre-Test Estimator, Mundlak Estimator, Fixed Effects Vector Decomposition

This paper proposes a new pre-test estimator of panel data models including time invariant variables based upon the Mundlak-Krishnakumar estimator and an "unrestricted" Hausman-Taylor estimator. The paper evaluates the biases of currently used restricted estimators, omitting the average-over-time of at least one endogenous time-varying explanatory variable. Repeated Between,

Introduction

We observe in sample time-invariant variables such as variables of origin, rare events or rarely changing variables which have been measured only once within the period of interest. Let us give a few examples of time invariant variables: colonial, legal or political system, international conflicts, institutional and governance indicators, initial gross domestic product per head when testing the convergence of incomes in growth regressions; geographical position for cross country data in gravity models of foreign trade and foreign direct investments; years of schooling, gender and race when testing wage income using survey data. These variables are often highly relevant in a theoretical model predicting correlations with a cross-sections and time-varying variable of interest. Because publishing results that do not reject the null hypothesis of no effect of these time invariant variable on the time varying variable of interest is generally not accepted by journal editors, applied researchers are likely to select estimators that reject the null hypothesis among the pool of available estimators.

We are then interested in drawing inference with respect to the statistical significance of these time-invariant variables (observed for N individuals) for explaining the variance of a time and individual varying variable (observed for N individuals during T periods). Since a time-invariant variable has no variance in the time direction, it can only explain the variance of a time and individual varying variable in its individual direction. The reason why it matters for inference is that "the effect of a random component can only be averaged out if the sample increases in the direction of that random (time or individual) component" [START_REF] Kelejian | Inference in Random Coefficient Panel Data Models: a Correction and Clarification of the Literature[END_REF], see also Hsiao 2003, pp. 51-53.). For example, the pooled ordinary least square (OLS) estimator with time-invariant explanatory variables leads to inference on time-invariant variables using the time dimension NT . But "we have to deal with the unfortunate fact that there is not quite so much information in N individuals observed T times as there is with NT individuals." (Johnston and Di Nardo (1997), pp. 395.). This leads to the suggestion of a number of possible estimators for the parameter values with different characteristics. The available estimators of time invariant variables using panel data or time series cross sections are found in [START_REF] Hausman | Panel data and unobservable individual effects[END_REF], [START_REF] Hsiao | Analysis of Panel Data[END_REF], [START_REF] Baltagi | Fixed Effects, Random Effects or Hausman-Taylor? A Pre-test Estimator[END_REF], [START_REF] Oaxaca | Fixed effects models with time invariant variables: a theoretical note[END_REF], [START_REF] Krishnakumar | Time Invariant Variables and Panel Data Models: A Generalised Frisch-Waugh Theorem and its Implications[END_REF] extension of [START_REF] Mundlak | On the pooling of time series and cross section data[END_REF] and Plümper and Troeger (2007) Fixed Effect Vector Decomposition Estimator (FEVD). For example, the FEVD estimator has been widely used over the recent years (for example, [START_REF] Blaydes | Rewarding Impatience: A response to Goodrich[END_REF] and [START_REF] Goodrich | A Comment on "Rewarding Impatience[END_REF]). By contrast, the Mundlak-Krishnakumar (2006) estimator remained unnoticed by practitionners.

This paper suggests that inference on time-invariant variables without a pre-test Mundlak-Krishnakumar estimator may lead to conclude wrongly that time-invariant variables are statistically significant. This paper includes four contributions:

(1) It proposes a pre-test estimator based upon the Mundlak-Krishnakumar estimator and a modified Hausman-Taylor estimator, extending the pre-test estimator proposed [START_REF] Baltagi | Fixed Effects, Random Effects or Hausman-Taylor? A Pre-test Estimator[END_REF].

(2) The Mundlak-Krishnakumar model is not currently used for doing inference on time invariant variables. It means that some currently published papers where the data generating process could be approximated by the Mundlak-Krishnakumar theoretical model report, for example, restricted estimators omitting the average over time of endogenous time-varying variables. The paper computes the determinants of the bias of the estimated parameters and of the estimated standard errors of using other estimators instead of the pre-test estimator. It presents an illustration of these omitted variable biases on a returns to schooling classic database [START_REF] Baltagi | On efficient estimation with panel data: an empirical comparison of instrumental variables estimators[END_REF], Cornwell and Rupert (1988), [START_REF] Baltagi | Econometric Analysis of Panel Data[END_REF]2009), [START_REF] Cameron | Microeconometrics Using Stata[END_REF], chapter 8).

(3 and 4) The paper explains why the [START_REF] Oaxaca | Fixed effects models with time invariant variables: a theoretical note[END_REF] and FEVD estimators may provide misleading estimated standard errors of the parameters of time-invariant variables. [START_REF] Mitze | Endogeneity in Panel Data Models with Time-Varying and Time-Fixed Regressors: To IV or not IV?[END_REF] and Breusch, Ward, Nguyen and Kompas (2010) found that FEVD estimated standard errors of the estimated parameters of time invariant variables were (too) small. The explanation that we propose is that the FEVD estimator uses the within root mean square error for the estimator of these standard errors. In so doing, it contradicts the Frisch-Waugh-Lovell theorem, where the same orthogonal projection matrix has to be used in the parameter estimate and in the estimator of its variance (following a similar argument, Greene (2010) refers to Aitken's theorem).

The overall picture is that when one uses this pre-test estimator, results may differ a lot from currently published inference which neglects a pre-test stage. There are two reasons for that. First, an omitted variable bias of the average over time of endogenous time-varying variables may occur. Secondly, the estimator of the standard error of the parameters of time invariant variables may be biased because of an excessive weight on within mean square error and on within degrees of freedom (depending on the time dimension T ). As the extent of these potential biases are not known when the data sets are not available for replication, this casts serious doubts on the emphasis put on inference on time invariant variable in panel data formerly published in academic journals.

This paper proceeds as follows. Section 2 defines our pre-test estimator based upon the Mundlak-Krishnakumar and the Hausman and Taylor estimator. Section 3 compares the outcome of different estimators using a typical time series cross section data set of gasoline demand. Section 4 analyses the differences of our pre-test estimator with the currently available estimators, in particular, the Oaxaca-Geisler, the two stage, the FEVD and the usual random effects estimators. Section 5 concludes.

A Pre-Test Estimator including Time-invariant Variables and Correlated Individual Effects

Mundlak-Krishnakumar Estimator

The model of time-series cross-sections regression estimates the following equation

y it = X it β + Z i γ + α i + ε it (2.1)
where y it denotes the endogenous variable, X it is a NT × k matrix of cross-sections time-series data, Z i is a NT × g matrix of time invariant variables, β and γ are k and g vectors of coefficients associated with time-varying and time-invariant observable variables respectively. Subscripts indicate variation over individuals (i = 1, ..., N)

and time (t = 1, ..., T ). Observations are ordered first by individual and then by time, so that α i and each column of Z i are NT vectors having blocks of T identical entries within each individual i = 1, ..., N. The disturbance ε it is assumed to be uncorrelated with the columns of (X it , Z i ,α i ) and has zero mean and constant variance σ 2 ε conditional on X it and Z i . The individual effect α i is assumed to be a time-invariant random variable, distributed independently across individuals with variance σ 2 α . The primary focus of the literature is the potential correlation of α i with the columns of (X it , Z i ).

Mundlak [1978] and [START_REF] Krishnakumar | Time Invariant Variables and Panel Data Models: A Generalised Frisch-Waugh Theorem and its Implications[END_REF] introduce an auxiliary regression which takes explicitly into account a linear relation between the explanatory variables and the individual random effects:

α i = X i. π + Z i φ + α M i (2.2)
where it is assumed that the disturbance α M i ∼ (0, σ 2 α M ) and where X i. is the average over time for each individual of time varying variables, such that

E ³ α M i | X i. ´= 0,
y it = X it β + Z i (γ + φ) + X i. π + α M i + ε it (2.3)
When E (α i | Z i ) = φ 6 = 0, one may only estimate the sum γ +φ and may not identify γ and φ separately, without additional prior information. The prior information the [START_REF] Hausman | Panel data and unobservable individual effects[END_REF] procedure uses is the ability to distinguish columns of X it and Z i which are asymptotically uncorrelated with α i from those which are not.

For fixed T , let

plim N→∞ 1 N X 0 1it α i = π 1 = 0, plim N→∞ 1 N X 0 2it α i = π 2 6 = 0 plim N→∞ 1 N Z 0 1i α i = φ 1 = 0, plim N→∞ 1 N Z 0 2i α i = φ 2 6 = 0 where X it = [X 1it , X 2it ] and Z i = [Z 1i , Z 2i ] are split into two sets of variables such that X 1 is NT × k 1 , X 2 is NT × k 2 , Z 1 is NT × g 1 , Z 2 is NT × g 2 with k 1 + k 2 = k
and g 1 + g 2 = g. X 1 and Z 1 are assumed exogenous and not correlated with α i and ε it , while X 2 and Z 2 are endogenous due to their correlation with α i but not with ε it . The pre-test estimator proposed in the next section differs from the Hausman and

Taylor [1981] estimator in that no prior information is assumed about how to split π, but that the information is obtained using t-tests on each estimated parameter b π of equations (2.2) or equation (2.3).

It will prove helpful to recall the menu of conventional estimators for (β, γ) in equation (2.3). Letting i T denote a T vector of ones, two orthogonal projection operators can be defined as: ) is always equal to zero: cov (y ity i. , x i. ) = 0 and cov (y ity i. , y i. ) = 0 (2.4)

B = I N ⊗ 1 T i T i 0 T , W = I NT -B
Hence, the total sum of squares SST of y it is the sum of total sum of squares SST W for the within transformed variable and T times the total sum of square SST B for the between transformed variable over N observations:

SST (y it ) = SST (y it -y i. ) + SST (y i. ) = SST W + T • SST B . (2.5)
Transform model (3) by the Within projection operator, we obtain:

y it -y i. = WX it β + WZ i (γ + φ) + WX i. π + Wα M i + Wε it (2.6) y it -y i. = (X it -X i. ) β + ε it -ε i. (2.7)
Least square estimates of β in the within transformed equation are Gauss-Markov for the transformed equation and define the within-groups estimator:

b β W = (X 0 it FX it ) -1 X 0 it Fy it with F = W. (2.8)
Since the columns of WX it are uncorrelated with Wε it , b β W is unbiased and consistent for β regardless of possible correlation between α i and the columns of X it or Z i . The OLS analysis of variance in the within subspace of observations (SSM is the sum of squares of the model, SSE is the sum of squares of the error) is:

SST W = SSM W + SSE W .
The sum of squared residuals (denoted form this equation can be used to obtain an unbiased and consistent estimate of σ 2 ε :

MSE W = c σ 2 ε = SSE W NT -N -k (2.9)
To make use of between-group variation, transform model (3) by the Between projection operator obtaining

y i. = X i. (β + π) + Z i (γ + φ) + α M i + ε i. When E (α i | X it , Z i ) = 0 = π = φ,
SST B = SSM B + SSE B .
The sum of squared residuals SSE B provides a biased and inconsistent estimator for

var (α i + ε i. ) = σ 2 α i + (1/T ) σ 2 ε when E (α i | X it , Z i ) 6
= 0 and unbiased and consistent when E (α i | X it , Z i ) = 0. Whenever we have consistent estimators for both β and γ, a consistent estimator for σ 2 α can be obtained with:

MSE B = c σ 2 α + c σ 2 ε T = SSE B N -k -g -1 .
(2.10)

Estimators that correspond to weighted average of the within and the between estimators are obtained as follows. Summing the within equation and the between equation multiplied by a parameter θ and repeated T times, which amounts to transform equation (3) applying the linear operator W+θB, we obtain:

y it -y i. + θy i. = (W+θB) X it β + θZ i. (γ + φ) + θX i. π + θα M i + (W+θB) ε it (2.11)
When π 6 = 0, that is when E (α i | X it ) 6 = 0, and when [START_REF] Mundlak | On the pooling of time series and cross section data[END_REF] and [START_REF] Krishnakumar | Time Invariant Variables and Panel Data Models: A Generalised Frisch-Waugh Theorem and its Implications[END_REF] prove that the best linear unbiased estimator is obtained for a value of θ corresponding to the generalized least square model (GLS) or random effect estimator, which leads to the following estimates:

E (α i | Z i ) = 0 = φ,
b β GLS = b β W , b π GLS = b β B -b β W and b γ GLS = b γ B , (2.12
)

var (b π GLS ) = var ³ b β B ´+ var ³ b β W ´and var (b γ GLS ) = var (b γ B ) (2.13)
A demonstration for the equalities on estimated variance-covariance matrix of estimated parameters is also found in Baltagi (2009, p.152-154). More precisely, [START_REF] Mundlak | On the pooling of time series and cross section data[END_REF] shows that (i) if α i is correlated with every column of X i.

σ b γ GLS = b σ b γ B = s SSE B N -k -g -1 1 q CSS(z j ) q 1 -R 2 A (z j ) , for j = 1, ..., g (2.14)
where CSS is the sum of squares corrected for the mean and 1 -R 2 A (z j ) is the "tolerance" (or the inverse of the variance inflation factor) and where R 2 A (z j ) is the coefficient of determination of the following auxiliary regression where the explanatory variable z i is explained by all other explanatory variables of regression (3) in the between dimension [START_REF] Stewart | Collinearity and Least Squares Regression[END_REF]):

z j = X i. π 0 + Z i,-j γ 0 + η 0 j
where η 0 j are disturbances, where π 0 and γ 0 are coefficients, and where Z i,-j is the matrix of time-invariant variables excluding the column of observations related to the variable z j . R 2 A (z i ) measures the effect of other explanatory variables on the estimated standard error of the estimated parameter of a given variable z i . When

R 2 A (z i
) is close to unity, there is a potential problem of near-multicollinearity. The root mean squared error of the between regression (including the degrees of freedom in the denominator) is

q SSE B N -k-g-1 .
The convergence of the estimator is obtained by increasing N to infinity, for T fixed. An increase of T does not lead the between estimator to converge. Inference on time-invariant explanatory variables do not depend on the number of observations in the time dimension T .

A Pre-test Estimator

The pre-test estimator consists of a first step using the Mundlak-Krishnakumar estimator (equation (2.3) using GLS method) and then perform specification tests of null

hypothesis H 0 : E (α i | X mit ) = π m = 0 against the alternative H 0 : E (α i | X mit ) = π m 6 = 0 for each of the time-varying explanatory variables 0 ≤ m ≤ k. If H 0 is rejected
for exactly k 2 variables, these variables are included in the subset X 2it assumed to be endogenous due to their correlation with α i but not with ε it , while the remaining k 1 variables are included in the subset X 1it . Hence, no prior information on how to split X it is required for this pre-test estimator. In practice, researchers priors on how to split X it and Z i based on theory may be misleading because of unobservable and omitted time-invariant variables, measurement errors, and so on.

The second step consists of using an unrestricted Hausman Taylor (1981) estimator using X 1it as instruments for Z 2 endogenous time invariant variables, but keeping the time invariant variables X 2i. . which are not used as instruments, in order to correct for their endogeneity (this is not usually done with the Hausman and Taylor estimator).

More precisely, if the specification tests rejects H 0 for k 2 ≥ 1 variables, it may still be possible to obtain consistent estimates of both β and γ in a second stage. Let

b d i = y i. -X i. b β W = ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´yit
be the NT vector of group means estimated from the within-groups residuals. Ex-panding this expression using equation 2.3 including only X 2i. leads to:

b d i = Z 1i. γ 1 + Z 2i. (γ 2 + φ 2 ) + X 2i. π 2 + α M i + ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´εit . (2.15)
Treating the last two terms as an unobservable mean zero disturbance, consider estimating γ from the above equation using N observations. If α i is correlated with the

columns of Z 2i , E (α i | Z 2i ) = φ 2 6
= 0, according to prior information, both OLS and GLS will be inconsistent estimates for γ. Consistent estimation is possible, however, if the columns of X 1it , uncorrelated with α i according to the non rejection of the null hypothesis of preliminary tests, provide sufficient instruments for the columns of Z i.

in equation (4.1). A necessary condition for identification of γ 2 and φ 2 is that k 1 ≥ g 2 :

there are at least as many exogenous time-varying variables as there are endogenous time-invariant variables (proposition 3.2 in Hausman and Taylor (1981)). When the condition k 1 ≥ g 2 is fulfilled, one may proceed to a second step for estimating γ knowing that the first step provided efficient estimates of β.

The two stage least squares (2SLS) estimator for γ in equation (4.1) is:

b γ II = ([Z 0 i , X 2i. ] P A [Z 0 i , X 2i. ]) -1 [Z 0 i , X 2i. ] 0 P A b d i (2.

16)

where A = [X 1it , Z 1i ] and P A is the orthogonal projection operator onto its column space. The sampling error is given by

b γ II -γ = ([Z 0 i , X 2i. ] P A [Z 0 i , X 2i. ]) -1 [Z 0 i , X 2i. ] 0 P A ³ α M i + ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´εit ánd
under the usual assumptions governing X it and Z i , the 2SLS estimator is consistent for γ, since for fixed T , plim N→∞

1 N A 0 α i = 0 and plim N→∞ 1 N X 0 it ε i = 0.
Having consistent estimates of β and, under the condition k 1 ≥ g 2 , γ, we can construct consistent estimators for the variance components. A consistent estimate of σ 2 ε can be derived from the within-group residuals in the first step c σ 2 ε = MSE W . Whenever we have consistent estimators for both β and γ, a consistent estimator for σ 2 α can be obtained.

Let

s 2 = (1/N ) ³ Y i. -X i. b β W -Z i b γ II -X 2i. b π 2,II ´0 ³ Y i. -X i. b β W -Z i b γ II -X 2i. b π 2,II
´:

then plim N→∞ s 2 =plim N→∞ 1 N (α i + ε i ) 0 (α i + ε i ) = σ 2 α + 1 T σ 2 ε so that s 2 a = s 2 -(1/T ) s 2 ε is consistent for s 2 a .
In the particular case when one assumes that α i is uncorrelated with the columns of

Z i : E (α i | Z i ) = 0 = φ,
and when all time varying variables reject the null hypothesis H 0 (k 1 = 0) using N observations to estimate the above equation in a second stage II, both OLS and GLS will be consistent estimates of γ identical to the between estimates

( b γ II,OLS = b γ B , b σ b γ II,OLS = b σ b γ B ) and with b π II,OLS = b β B -b β W . Baltagi,
X it = [X 1 , X 2 ], (3) it includes X 2i. π 2 which
may change widely the estimates of γ in the second stage, (4) when 0 ≤ k 1 < g 2 (not enough exogenous X 1it to instrument endogenous Z 2i ) or respectively when g 2 = 0

and k 1 = k (no endogenous Z 2i , and all time varying variables are endogenous X 1it ) it reports a biased (respectively unbiased) estimators for γ, instead of no estimates at all.

Because classical statistical hypothesis testing implies non-zero probabilities of type I error (p-value) and type II error (one minus the power of the test), a pre-test estimator cannot perform as well as an estimator where the researcher exactly knows the "true" model based on "true prior information" before testing. This only occurs when the number of observations tends to infinity, so that the probabilities of both types of errors tends to zero. [START_REF] Mundlak | On the pooling of time series and cross section data[END_REF] considers the following estimators for his theoretical model: pooled OLS (θ = 1), between, within (θ = 0), Generalized least square (GLS) or random effects estimator (for an estimated value of θ), restricted generalized least square with the restriction π = 0 (that is applying the usual random effect model whereas the true model is such that: π 6 = 0), restricted pooled OLS with the restriction π = 0. The reasons for considering restricted estimators are twofold. Firstly, these restrictions are likely to decrease the variance of the estimators although these restricted estimators are generally biased. There is therefore a trade off between bias and variance and the choice of an estimator depends on the weights to be assigned to the two components.

Secondly, the Krishnakumar-Mundlak is not currently used for doing inference on time invariant variables. It means that some currently published papers where the data generating process could be approximated by the Mundlak theoretical model report Mundlak's restricted estimators. It is interesting to evaluate the bias of the estimated parameters and of the estimated standard errors of the time invariant variables for the restricted models. We therefore do a similar investigation as Mundlak (1978) did on retricted estimators for time varying variables, and we focus on the estimators related to time-invariant variables when our pre-test is estimator is consistent. In tables 1 and 2, the dependent variable log(wage) is explained by nine timevarying variables (X it 's) and three time-invariant variables (Z i 's). This lack of precision is partly due to the fact that X * 1 consists of three weak instruments explaining altogether only 12% of the variance of years of education. The pre-test estimator leads to choose highly exogenous instruments X 1 . If the variables Z 2 are strongly correlated with unobservable heterogeneity α i , then the exogenous instruments X 1 will be weakly correlated with Z 2 . Because the pre-test estimator leads to choose weak instruments, the Mundlak estimator is likely to be as reliable as the Hausman-Taylor estimator for this data set.

A Return to Schooling Illustration of the Unrestricted Model

σ b γ GLS = b σ b γ B = q SSE B N-k-g-1 q CSS(z i ) q 1 -R 2 A (z j ) = q 42.

Repeated Between Estimator

We present this estimator because repeating T times the time invariant observations turns out to be the major component of the increase of the t-statistic of time-invariant variables using panel data for several other estimators. However, [START_REF] Kelejian | Inference in Random Coefficient Panel Data Models: a Correction and Clarification of the Literature[END_REF] argue "the effect of a random component can only be averaged out if the sample increases in the direction of that random (time or individual) component". By contrast, Oaxaca and Gleiser (2003) assume that the consistency of the estimator of a parameter of a time-invariant explanatory variable "depends on the time series observations approaching infinity". The estimated parameters of time-invariant variables are the same with the repeated between and the between estimator (γ B = γ RB ). The between estimator of the standard error with observations repeated T times is equal to (the subscript for this estimator is RP B):

b σ RB b γ = q T •SSE B NT -k-g-1 q T • CSS(z j ) 1 q 1 -R 2 A (z j ) = s N -k -g -1 NT -k -g -1 b σ B b γ = b σ B b γ

.

Inference uses NT -k -g -1 degrees of freedom instead of N -k -g -1 in the between estimator. The coefficient of determination of the auxiliary regression R 2 A (z j ) does not change in the between or repeated between samples. The estimated standard error of the estimated parameter of the time invariant variables is divided by q 4165-13 595-13 = 2.67 in table 1. As the parameter estimate is the same as γ B , the repeated between t RBstatistic amounts to multiply the between t B statistics by the following factor:

t RB = c γ B b σ RB c γB = s NT -k -g -1 N -k -g -1 c γ B b σ B b γ = s NT -k -g -1 N -k -g -1 t B
When N is large, the t-statistic of the repeated between model is multiplied by around √ T (say by 2 when T = 4 and by 5 when T = 25) with respect to the between model.

Pooled OLS Unrestricted Estimator

The Pooled OLS (θ = 1) estimator ignores the random effects. It leads to the same parameter estimates as the Mundlak estimator, but not the same standard error estimates. The reason why we analyse the Pooled OLS estimator for time-invariant explanatory variables is that it is still used by some applied researchers. When test-ing does not reject the existence of random effects and when the Hausman test leads to reject the random effect models with respect to the "correlated random effects" model, some applied researchers present Pooled OLS estimates with time invariant explanatory variables, on the ground that the required within estimates wipes out time-invariant variables. For example, [START_REF] Oaxaca | Fixed effects models with time invariant variables: a theoretical note[END_REF] evaluate the consistency of the OLS estimator of a parameter of a time-invariant explanatory variable.

The OLS estimator uses the irrelevant NTkg -1 degrees of freedom for time-invariant variables. The estimated standard error is larger than for the repeated between estimator, because the RMSE of the OLS estimator is larger since it includes the sum of squares of the error of the within model:

SSE OLS = T • SSE B + SSE W (the last equality presents table 1 values): b σ OLS b γ = q T •SSE B +SSE W NT -k-g-1 q T • CSS(z i ) 1 q 1 -R 2 A (z j ) = s 1 + SSE W T • SSE B s N -k -g -1 NT -k -g -1 b σ B b γ . = s 1 + 84.11480 7 • 42.07257 = 1.1338 • b σ B b γ 2.671 = b σ B b γ 2.33 . R 2 A (z j )
y it -y i. + b θy i. with b θ = RMSE W √ T • RMSE B = 0.14071 √ 7 • 0.26887 = 0.1978.
The between regression includes the same variables (the average-over-time of all explanatory variables. Hence, the root mean squared error of the between estimator RMSE B is the same in both models. The restricted random effect model faces an omitted variable bias when one rejects the null hypothesis b β B = b β W , for at least one time-varying variable. In the wages example, the bias with the random effect model is (γ GLS denotes the Mundlak GLS estimator whereas γ RGLS denotes the Mundlak restricted GLS estimator):

b γ RGLS = b γ B + j=k X j=1 ³ b β B -b β W ´b β b θx i. / b
θz ii = 0.06966 + 0.03044 = 0.10 with β x i. /z i given by the auxiliary regression using OLS on quasi demeaned variables (table 2, column 1):

b θx i. = β x i. /z i b θz i + β x i. /x it ³ b θx i. + x it -x i. ´+ b θβ 0 + ε it .
The omitted variable bias on the parameters of the time-invariant variables is large for ED (the return of education parameter nearly doubled: it is multiplied by 1.93β B ) and it is smaller for BLK (1.33β B )and FEM (1.06β B ). The omitted variable bias on the t-statistic is 1.87t B for ED, 1.03t B for BLK, 1.14t B for FEM.

Restricted Hausman-Taylor Estimator

In table 2, the column HT(X * * 1 ) presents Baltagi and Khanti-Akom (1990) estimates using the standard Hausman Taylor procedure, with the set of instruments X * * 1 and omitting the average-over-time of the endogenous time-varying variables X * * 2 . The return to schooling is 0.137 with a confidence interval [0.116, 0.158].

The rise of the return to schooling with respect to the between estimate is only partly due to "instrumental variables" correcting the endogeneity of schooling with the individual effect. Around half of the rise is given by the omission of time-invariant variables which corrects the endogeneity of time varying regressors with the individual effect. The other component of the rise of the parameter is related to the strong instrument "blue collar occupation" which is strongly correlated with the number of years of education (its correlation coefficient around 0.45). This is likely to drive the results emphasized by [START_REF] Boumahdi | Instrument relevance and efficient estimation with panel data[END_REF] for their measure of instruments relevance on this data set. They find that the gain in efficiency across the various sets of instruments seems to be in the education variable. However, "blue collar occupation" is endogenous, with a relatively large difference between its estimated

parameters b β B -b β W with
d i = y i. -X i. b β W = ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´yit
be the N vector of group means estimated from the within-groups residuals. Expanding this expression leads to:

b d i = Z 1i. γ 1 + Z 2i. (γ 2 + φ 2 ) + α M i + ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´εit + X i. ³ b β B -b β W (4.1)
Treating the last two terms as an unobservable mean zero disturbance, consider estimating γ from the above equation using N observations, without taking into account that the disturbances of the restricted between includes this omitted term:

X i. ³ b β B -b β W ´.
There is then an omitted variable bias on the parameter of the time-

invariant variable: b γ II-RB = b γ B + j=k X j=1 ³ b β B -b β W ´b β x j. /z i
with β x i. /z i estimated using the following auxiliary regressions in the between dimension:

x i. = β x i. /z i z i + β 0,x i. /z i + ε i,x i. /z i .
In this two stage restricted between estimator, the omitted variable bias on the estimated standard error of the estimated parameter of the time-invariant variables leads to differences with respect to the between estimator:

b σ 2 b γ II-RB = s SSE II-RB N -g -1 1 q CSS(z i ) q 1 -R 2 A,RB (z j ) 6 = b σ b γ B = s SSE B N -k -g -1 1 q CSS(z i ) q 1 -R 2 A (z j )
First, the sum of squares of errors SSE II-RB is larger than the SSE B because the model constrains the parameters of the time-varying variables to their within estimate which may not minimize the between sum of squares of errors. However, this increase of the estimated standard error may be offset for two reasons due to the fact the averages of k time-varying variables are on the left hand side of the equation. First, the degrees of freedom increase by k. This decreases the root mean squared error.

Second, the variance inflation factor decreases because

R 2 A,RB (z j ) < R 2 A (z j ). R 2 A (z j )
is the coefficient of determination of an auxiliary regression where the time-invariant variable is correlated with the other g -1 time-invariant explanatory variables on the right hand side of the equation. In the between estimator, R 2 A (z j ) is the coefficient of determination of an auxiliary regression where the time-invariant variable is correlated with the other g-1 time-invariant explanatory variables and k averages of time-varying explanatory variables. As the number of explanatory variables increases by k, one has

R 2 A,RB (z j ) < R 2 A (z j ).
With respect to the estimated parameter of the non restricted between, the estimated parameter of the restricted between are multiplied by 0.89 for BLK, by 1.23

for ED and by 1.42 for FEM. With respect to the t-statistics of the non restricted between, the t-statistics of the restricted between are multiplied by 0.77 for BLK, by 1.47 for ED and by 1.86 for FEM. [START_REF] Oaxaca | Fixed effects models with time invariant variables: a theoretical note[END_REF] propose an alternative correction of the omitted variable bias of the two-stage restricted between approach than the Mundlak-Krishnakumar estimator. They estimate the repeated restricted between with a two stage generalized least square estimator using a covariance matrix that takes into account that the disturbances of the restricted between includes this omitted term:

X i. ³ b β B -b β W
´.

The first drawback of the [START_REF] Oaxaca | Fixed effects models with time invariant variables: a theoretical note[END_REF] estimator is that the second step between should not be estimated using observations repeated T times. The second drawback is that it is simpler to include explicitly the variables X i. in the Mundlak-Krishnakumar estimator, than to compute a two stages GLS estimator with a specific covariance matrix.

Three-stage FEVD Unrestricted and Restricted Estimator

The FEVD estimator adds a third stage to the previous two stage estimator. It is assumed that time-invariant variables are not correlated with the random individual effects: E (α i | Z i ) = φ = 0 and three sets of alternative assumptions can be dealt which is consistent with the assumption E (α i | X it ) = π = 0. Note that if α i is random instead of being fixed, the "usual" GLS estimator is consistent, which is not the case of the first stage "within-groups" estimator for the FEVD.

with (E (α i | X it ) = π = 0, or E (α i | X 2it ) = π 2 6 = 0 or E (α i | X it ) = π 6 = 0 (
Let us denote b ε i,B the between residuals of stage II restricted between:

b ε i,B = b α M i + ³ B -X i. (X 0 it WX it ) -1 X 0 it W ´b ε it +X i. ³ b β B -b β W ´= y i. -X i. b β W -Z i. b γ II-RB
The third stage of the FEVD estimator is an OLS regression which includes the residual of the second stage repeated between regression denoted b ε i,B with a parameter δ estimated or constrained to one.

y it = X it β + Z i γ + X i. π + b ε i,B • δ + η it,III .
Expanding this expression yields:

y it -y i. = Wy it = (X it -X i. ) β -X i. ³ β -δ b β W ´+ Z i (γ -δb γ II-RB ) + η it,III (4.2)
Wy it is orthogonal in the sample to all time-invariant variables, so that OLS estimates This analysis explains why Kristensen and Wawro 2007 (footnote p.22) found that the FEVD estimated standard errors were relatively "too small" using Monte Carlo simulations with respect to other estimators. The FEVD estimator of the standard error is (the last equality refers to the unrestricted FEVD case, table 1

are: b γ = b δb γ B , b π = b δ ³ b β B -b β W ´, b β = b δβ W , b β = b β W , b δ = 1,
column 6): b σ b γ F EV D = q SSE W NT -N-k q T • CSS(z i ) 1 q 1 -R 2 A (z j ) = s N -k -g -1 NT -N -k s SSE W T • SSE B b σ B b γ = s 595 -13 4165 -595 -9 s 84.11480 7 • 42.07257 b σ B b γ = 1 2.4736 • 1 1.8712 • b σ B b γ = b σ B b γ 4.6286
The FEVD estimated standard error of a time-invariant variable is usually biased downwards for two reasons:

-It uses NT -Nk degrees of freedom (with N the number of individuals, k the number of time-varying explanatory variables, with T the number of periods) instead of Nkg -1 degrees of freedom.

-It multiplies the repeated between estimator of the standard error by a positive factor q SSE W T •SSE B which can be much smaller than one when T increases. Taylor models with perfect and imperfect knowledge about the underlying correlation between right hand side variables and residual term. Nonetheless, concluding that FEVD is "more efficient" is misleading, because the computation of the estimated standard errors of estimated parameters contradicts Frisch-Waugh-Lovell theorem.

Conclusion

Our theoretical and empirical investigation on inference on time-invariant variables shows that a pre-test estimator based upon the Mundlak-Krisnakumar and a modified

Hausman-Taylor estimator should be used. Furthermore, the procedures are already programmed and available in all econometric softwares. The first stage consists of the usual random effects GLS estimator including all the variables X i. . The second stage uses instrumental variables X 1i. keeping the X 2i. as explanatory variables, using instrumental variables estimator or the Hausman-Taylor estimator procedures.

An example shows first that a time-invariant variable is not statistically significant for some estimators and highly statistically significant with other estimators.

π

  and φ are k and g vectors of coefficients associated with the average over time of time-varying variables and time-invariant observable variables, respectively. Clearly π = 0 and φ = 0 if and only if the time varying and time invariant variables are uncorrelated with the random individual effects. Combining the auxiliary regression with the initial regression yields:

  least squares estimates of β and γ in the between transformed regression define the between groups estimators estimators (denoted b β B and b γ B ) which are unbiased and consistent for β and γ, using N observa-tions. When E (α i | X it ) 6 = 0 or when π 6 = 0 , b β B is biased and inconsistent. When E (α i | Z i ) 6 = 0 or when φ 6 = 0, b γ B is biased and inconsistent. The OLS analysis of variance in the between subspace of observations is:

(π 6

 6 = 0, with k 1 = 0 and k 2 = k), the Gauss-Markov estimator for β is the within groups estimator b β W , and (ii) if α i is uncorrelated with every column of X i. (π = 0 in the true model, with k 1 = k and k 2 = 0), the Gauss-Markov estimator for β is the "usual" GLS estimator b β GLS , (iii) Using the "usual" GLS estimate assuming the restriction π = 0 when the true model is such that π 6 = 0 is denoted "'restricted GLS" (RGLS) and leads to biased estimate b β RGLS 6 = b β GLS = b β W , because of the k omitted variables X i. bias. The estimator b π GLS is equal to the difference between the between estimator and the within estimator b β 1,B -b β 1,W . For the time-invariant variables, the estimator b γ GLS is exactly the between estimator b γ B . The standard error estimator, b σ b γ GLS , of the estimated parameter using the Mundlak GLS model is exactly the same as the one of the between regression (see appendix 1). This estimator -even though the data includes NT observations -takes into account that for time-invariant explanatory variables only N observations should be used. It has Nkg -1 degrees of freedom: b

3. 1 .

 1 Pre-test Estimator The pre-test estimator leads to dramatic changes with respect to alternative estimators of time invariant variables in panel data. This is demonstrated empirically for a return to schooling example based on a panel of 595 individuals observed over the period 1976-1982 drawn from the Panel Study of Income Dynamics (PSID) (Baltagi and Khanti-Akom (1990), Cornwell and Rupert (1988)). First, the tests based on Mundlak-Krishnakumar specification lead to alternative choices of the instruments than the ones chosen by Baltagi and Khanti-Akom (1990) and Cornwell and Rupert (1988). Second, when the average over time of endogenous time varying variables are not omitted in the Hausman-Taylor estimates, the estimates of the return to schooling are much lower than the ones found by Baltagi Khanti-Akom [1990] and Cornwell and Rupert [1988].

  is the average number of years of experience in the sample of years 1976-1982 for each individual i and t is a time trend (= -3 for 1976,...,= +3 for 1982). This individual time trend accounts for 65% of the within variance of the dependent variable log(wage), along with a correlation coefficient of 0.8. Once this variable is included in the model, the partial R2 contribution of the eight other within transformed variables is very small (the R2 increases by only 0.74%). The estimated coefficients of blue collar occupation (OCC), of marital status (MS), of living in the south of the U.S.A. (SOUTH), of weeks worked (WKS), of industry (IND), of the standard metropolitan statistical area (SMSA), of belonging to an union (UNION) and of EXP 2 are statistically significant. The between estimator includes estimates of the three time invariant variables: female workers (FEM), black worker (BLK), and years of schooling (ED). Degrees of freedom are 582. All the estimated coefficients are statistically significant. The 95% confidence interval for the return to school is [0.045, 0.057]. The standard error estimator, b σ b γ GLS , of the estimated parameter using the Mundlak-Krishnakumar GLS model is exactly the same as the one of the between regression (table 1, column 2): b

4. 3 .

 3 Two-stage Restricted Between A common practice consists of a two-stage restricted between (denoted II -RB) estimator of time-invariant variables, with the restriction b β = b β W for the averageover-time of time-varying variables in a between regression. Let b

  all time-varying variables are endogenous). Plümper and Troeger (2007) use a restricted FEVD estimator with a second stage and a third stage omitting variables X 2i. or X i.

  and b η it,III = b ε it -b ε i. . The residuals of the third step regression are exactly the within-groups regression residuals. The only change of the three stage FEVD estimator with respect to the two stage estimator is related to the estimated standard error of the parameters γ of time invariant variables Z i . The FEVD estimator of the standard error of the estimated parameter of a time-invariant variable amounts to substitute the mean squared error in the between dimension (MSE B = SSE B / (Nkg -1)) by the mean squared error in the within dimension (MSE W = SSE W / (NT -Nk)). The estimated parameters b γ II-RB are related to the projection in the between subspace of observations, but their estimated standard errors are related to the projection in the within subspace of observations, which is orthogonal to the between subspace. This contradicts the Frisch-Waugh-Lovell theorem: the same orthogonal projection matrix has to be used in the parameter estimate and in the estimator of its variance. The Frisch-Waugh-Lovell decomposition of the Between and the Within dimension holds with any of these assumptions: effects are fixed or random, or the timevarying variables are correlated or not with the random effect do not matter at all.As the residuals in the between dimension are excluded from the computation of the variance of the parameters, the potential correlations of the time-invariant variable with the individual random effect (which are in the between dimension) are excluded from the computation of the FEVD variance of the parameters. As well, adding a large number of time-invariant variables in the regression, in particular when they are near-multicollinear with the time-invariant variable of interest, changes the estimated parameter, but does not increase the FEVD estimated standard error of the estimated parameter of the time-invariant variable of interest. So the FEVD wipes out nearmulticollinearity problems in the time-invariant dimension, which is also "practical"for reducing estimated standard errors.

The Mundlak -

 Mundlak Krishnakumar regression reports within estimates and between estimates, with tests of the null hypothesis b β B -b β W = 0 for each time-varying explanatory variable. In the case where at least one (but not all) x i. have estimated parameter b β B -b β W small and not significantly different from zero, the Mundlak estimator suggests which of the x i. is exogenous. Then, these exogenous x i. are the ones to be used as instrumental variables in a variant of the Hausman and Taylor [1981] instrumental variables estimator, including the average over time of endogenous time-varying vari-ables. Using this pre-test estimator, one does not rely on subjective priors for deciding which time-varying variables are or are not endogenous in the Hausman and Taylor estimator. One may try further Chamberlain (1984) estimators where the correlation of explanatory variables with the random individual effects is not only contemporaneous, but can also be related to leads and lags of explanatory variables. It is not necessarily the failure of a model that the outcome of the tests is: ¯b β B ¯> b β W = 0. Some variables have zero or small correlation coefficients in the within subspace, and large correlation coefficients in the between dimension, just because it is the statistical information included in these rarely changing or time-invariant variables during the period of observations in available data sets. If we turn to the case of rarely changing variables instead of time-invariant variables, these variables are characterized by small within correlation coefficients with the explained variable and with all other explanatory variables. By contrast, they may have high between correlation coefficients with the explained variable and with all other explanatory variables. The variables are typically such that the null hypothesis b β B -b β W = 0 is rejected. The Mundlak estimator provides directly the t-test related to this null hypothesis. It is of particular interest for those variables, which are likely to be endogenous and significant variables in the Mundlak approach. If ever they do not have large between correlation coefficients, then the joint hypothesis b β B = b β W = 0 may not be rejected.

  By it = y i. . Similarly, W produces a vector of deviation from group means: i.e., Wy it = y ity i. . For any time-invariant vector of observations, Wz i = O. The orthogonality of the two projection operators is such BW = WB = O NT . The space

	which are idempotent matrices of rank N and NT -N respectively. With data grouped by individuals, B transforms a vector of observations into a vector of group means, observations transformed by the between operator (dimension N) and the subspace of observations transformed by the within operator (dimension NT -N). In other words, the covariance of the between transformed variable (time average for each individual i, with t a time index and T observations) x i. = 1 T P t=T t=1 x it with respect to a within such that of observations (R NT ) has an orthogonal decomposition in two spaces: the subspace of transformation of another variable (y it -y i.

  Bresson and Pirotte (BBP)(2003) pre-test estimator reverts to the usual random effects estimator if the standard Hausman test based on the "within-groups" versus the random effects estimators is not rejected. It reverts to the HT estimator if the choice of strictly exogenous regressors is not rejected by a second Hausman test based on the difference between the "within-groups" and HT estimators. Otherwise, this pre-test estimator reverts to the "within-groups" estimator.

	The pre-test estimator differs from BBP (2003) in that (1) it assumes Mundlak
	auxiliary regression, (2) it tests the endogeneity of each time-varying variables instead
	of relying on prior information for splitting

  Table 2 includes six estimates of the restricted model: step II restricted between, step II repeated restricted between, GLS, OLS, FEVD and HT(X * * 1 ), where the set of exogenous time varying variables X * * 1 is the one chosen by Baltagi and Khanti-Akom (1990). 1 For the unrestricted model, we observe the following. The estimates of the all parameters are identical for all estimators, except for the average over time of a timevarying variable in the Between or Repeated Between estimators which exclude the time-and individual-varying variable. For the time-and individual-varying variable, Mundlak-Krishnakumar, LSDV and FEVD estimators exhibit the same estimated standard deviation and value of the t-statistics, which differs from OLS. But, for the time-invariant variable, the estimated standard deviation and the values of the tstatistic vary dramatically. For the Between and the Mundlak-Krishnakumar estimators, the estimated standard deviations of time invariant variables are equal, whereas the estimated standard deviations for the repeated between estimate are 2.67 times smaller, the ones for the OLS estimate are 2.33 times smaller, and the ones for the FEVD estimate are 4.63 smaller. It seems, and this is confirmed later on by our theoretical analysis, that the standard deviation of the OLS and FEVD is systematically smaller and the value of the t-statistic systematically larger than for the Mundlak estimator. If this is the case, inference with the OLS and FEVD estimator will less often reject the null-hypothesis that γ is equal to zero, and more often conclude that the time-invariant variable has an effect on the endogenous variable statistically different from zero. In table 1, the within estimator, which is consistent, serves as a benchmark for the Hausman (1978) test. The experience variable (EXP) is a time trend which differs only in levels for each individuals, i.e. EXP it = EXP i. + t where EXP i.

[START_REF] Baltagi | On efficient estimation with panel data: an empirical comparison of instrumental variables estimators[END_REF]

and

Cornwell and Rupert (1988)

assumed that four of the X it variables and two of the Z i variables are uncorrelated with the individual effects. These are denoted by X 1 and Z 1 , respectively. They are listed at the bottom of table 2. The remaining X it and Z i are correlated with the individual effects. They are denoted by X 2 and Z 2 , respectively. Table

1

includes estimates of the unrestricted model using six estimators: within, between, Mundlak-Krishnakumar, pre-test unrestricted Hausman-Taylor, repeated between (RB), ordinary least square (OLS) and FEVD.

  is the determination coefficient of the auxiliary regression. It does not change with respect to the between, repeated between and OLS estimators. In table1example, as T = 7 is small, the OLS t OLS statistics (equal to 2.33 times the between t B statistics) is close to the repeated between t RB statistics (equal to 2.67 times the between t B statistics). However, as the time-invariant variables do not explain any variance of the within transformed explained variable, taking into account SSE W is

	irrelevant for doing inference on time-invariant variables.T
	4. The Omitted Variable Bias of the Restricted Models
	4.1. Restricted Random Effect Estimator
	Both, the restricted random effects and the Mundlak-Krishnakumar estimator use the
	same weight b θ for computing quasi demeaned variables (the second equality refers to
	our example) using the Swamy and Arora method (1972):

  a relatively large t statistics (t = 4.57) in the Mundlak equation. Unobservable ability affect primarily the number of years of education, which determine the outcome of a "blue collar" occupation.

Let us finally remark that Angrist and Pischke (2009, p.64-68) do not recommend to include the time varying variable OCC it (occupation) when investigating the causal link between education and wages, because it is an outcome variable which occurs at a later stage than education.

  The combination of the potential omitted variable bias of the parameter estimate of the time invariant variable ED in the restricted between (step II) and of the bias of the estimated standard error implies a t-statistics equal to t RF EV D = 5.69 • t B times the t-statistics of the Mundlak model for the restricted FEVD (table 2, column 6). For the unrestricted FEVD (table 1, column 6), only the bias of the estimated standard error matters, so that the t F EV D = 4.62 • t B times the t-statistics of the Mundlak model for ED. For the time invariant BLK, t RF EV D = 3.6 • t B and for FEM, t RF EV D = 8.6 • t B . Note that when T corresponds to an aggregate annual time series

	(T = 25,	√ T = 5), the increase of t-statistics is effect is likely to double with respect
	to the wages example (T = 7,	√ T = 2.6). Mitze [2009] also finds that the FEVD
	estimator tends to have a smaller root mean square error (rmse) than both Hausman-

Table 1 . Dependent variable: log wage: Unrestricted Models.

 1 FEM, BLK), X1*=(SOUTH, IND, EXP2) using pre-test selection.

			Within	Between GLS Mundlak Pre-Test	Repeated	OLS	FEVD
						HT(X1*)	Between
	Constant		-	5.12	5.12	5.18	5.12	5.12	-
				(0.203)	(0.203)	(0.861)	(0.076)	(0.08)
	EXP	(it) 0.113	-	0.113	0.114	-	0.113	0.113
			(0.003)		(0.003)	(.00226)		(0.005)	(0.003)
	EXP2	(it) -0.00042	-	-0.00042	-0.000449	-	-0.00042	-0.00042
			(0.00005)		(0.00005)	(0.0000485)		(0.00011)	(0.00005)
	WKS	(it) 0.00084	-	0.00084	0.000802	-	0.00084	0.00084
			(0.0006)		(0.0006)	(0.0006)		(0.0011)	(0.0006)
	SOUTH	(it) -0.0019	-	-0.0019	-0.0388	-	-0.0019	-0.0019
			(0.035)		(0.035)	(0.0272)		(0.0680)	(0.035)
	SMSA	(it) -0.0425	-	-0.0425	-0.0405	-	-0.0425	-0.0425
			(0.0194)		(0.0194)	(0.0193)		(0.0385)	(0.0194)
	MS	(it) -0.0297	-	-0.0297	-.0307	-	-0.0297	-0.0297
			(0.019)		(0.019)	(.0189)		(0.037)	(0.019)
	OCC	(it) -0.0215	-	-0.0215	-0.0222	-	-0.0215	-0.0215
			(0.014)		(0.014)	(0.0137)		(0.027)	(0.014)
	IND	(it) 0.0192	-	0.0192	0.0297	-	0.0192	0.0192
			(0.0154)		(0.0154)	(0.0144)		(0.031)	(0.0154)
	UNION	(it) 0.0328	-	0.0328	0.0325	-	0.0328	0.0328
			(0.0149)		(0.0149)	(0.0148)		(0.0296)	(0.0149)
	EXP(i)	(i.) -	0.0319	-0.0813	-0.0875	0.0319	-0.0813	-0.0813
				(0.0048)	(0.0054)	(.00272)	(0.0018)	(0.0053)	(0.0012)
	EXP2(i)	(i.) -	-0.00057	-0.00015	-	-0.00057	-0.00015	-0.00015
				(0.00011)	(0.00012)		(0.00004)	(0.00012)	(0.000026)
	WKS(i)	(i.) -	0.0092	0.00835	0.00874	0.0092	0.00835	0.00835
				(0.0036)	(0.0037)	(0.00448)	(0.0013)	(0.0019)	(0.0080)
	SOUTH(i) (i.) -	-0.057	-0.055	-	-0.057	-0.055	-0.055
				(0.026)	(0.043)		(0.0097)	(0.069)	(0.0091)
	SMSA(i) (i.) -	0.176	0.218	0.218	0.176	0.218	0.218
				(0.026)	(0.032)	(0.0386)	(0.0096)	(0.040)	(0.0069)
	MS(i)	(i.) -	0.115	0.145	0.151	0.115	0.145	0.145
				(0.048)	(0.051)	(0.0542)	(0.0168)	(0.042)	(0.011)
	OCC(i)	(i.) -	-0.168	-0.146	-0.152	-0.168	-0.146	-0.146
				(0.034)	(0.0365)	(0.169)	(0.0126)	(0.031)	(0.0079)
	IND(i)	(i.) -	0.058	0.0387	-	0.058	0.0387	0.0387
				(0.026)	(0.0298)		(0.0095)	(0.032)	(0.0064)
	UNION(i) (i.) -	0.109	0.0763	0.0835	0.109	0.0763	0.0763
				(0.029)	(0.0328)	(0.0381)	(0.0109)	(0.0321)	(0.0071)
	FEM(i)	(i)	-	-0.317	-0.317	-0.319	-0.317	-0.317	-0.317
				(0.055)	(0.055)	(0.0601)	(0.0205)	(0.023)	(0.012)
	BLK(i)	(i)	-	-0.158	-0.158	-0.165	-0.158	-0.158	-0.158
				(0.045)	(0.045)	(0.0499)	(0.0168)	(0.019)	(0.0097)
	ED(i)	(i)	-	0.0515	0.0515	0.0489	0.0515	0.0515	0.0515
				(0.00555	(0.00555)	(0.0469)	(0.00208)	(0.00235)	(0.0012)
	D.of F.		3561	582	3561 and 582	3561 and 582 4152	4143	3561
	RMSE		0.141	0.269	0.152		0.359	0.152	0.152
	R2		0.907	0.544	0.645		0.544	0.645	0.645
	Z1=(

Table 2 . Dependent variable: log wage: Restricted Models.

 2

			R-Between	R-GLS	R-HT (X1** Repeated	R-OLS	R-FEVD
			(step II)			R-Between
	Constant		5.92	4.264	2.913	5.92	5.25	-
			(0.061)	(0.098)	(0.283)	(0.023)	(0.07)
	EXP	(it)	-	0.082	0.113	-	0.040	0.113
				(0.003)	(0.019)		(0.002)	(0.003)
	EXP2	(it)	-	-0.0008	-0.000419	-	-0.0007	-0.00042
				(0.00006)	(0.000055)		(0.00005)	(0.00005)
	WKS	(it)	-	0.00084	0.00084	-	0.0042	0.00084
				(0.0008)	(0.0006)		(0.0011)	(0.0006)
	SOUTH	(it)	-	-0.0017	-0.0074	-	-0.0556	-0.0019
				(0.027)	(0.032)		(0.012)	(0.035)
	SMSA	(it)	-	-0.014	-0.0418	-	0.151	-0.0425
				(0.020)	(0.0189)		(0.120)	(0.0194)
	MS	(it)	-	-0.075	-0.0298	-	0.048	-0.0297
				(0.023)	(0.019)		(0.020)	(0.019)
	OCC	(it)	-	-0.050	-0.0207	-	-0.140	-0.0215
				(0.017)	(0.014)		(0.014)	(0.014)
	IND	(it)	-	0.004	0.0136	-	0.046	0.0192
				(0.017)	(0.0152)		(0.011)	(0.0154)
	UNION	(it)	-	0.063	0.0328	-	0.092	0.0328
				(0.017)	(0.0149)		(0.013)	(0.0149)
	FEM(i)	(i)	-0.449	-0.339	-0.131	-0.449	-0.367	-0.449
			(0.041)	(0.051)	(0.127)	(0.016)	(0.025)	(0.009)
	BLK(i)	(i)	-0.141	-0.210	-0.285	-0.141	-0.167	-0.141
			(0.051)	(0.058)	(0.155)	(0.019)	(0.022)	(0.011)
	ED(i)	(i)	0.0635	0.100	0.137	0.0635	0.0567	0.0635
			(0.0046)	(0.006)	(0.021)	(0.0017)	(0.0026)	(0.001)
	D.of F.		591	4152	3561 and 591 4161	4152	3561
	RMSE		0.313	0.1896		0.313	0.349
	R2		0.366	0.428		0.366	0.428
	Z1=(FEM, BLK), X1

**=(SOUTH, IND, OCC, SMSA), with t statistics for Mundlak estimation in table 1: OCC(i): t=-4.01, SMSA(i): t=6.77. R-OLS in Baltagi (2009), table 2.6, p.27. R-GLS in Baltagi (2009) table 2.8, p.28, R-HT (X1**) in Baltagi (2009), table 7.2, p.157.

Table 3 : t-statistics of time invariant variables, as a proportion of the t-statistics of the between estimator.

 3

			Between	GLS	Pre-Test	Repeated	OLS	FEVD
				Mundlak	HT(X1*)	Between
	FEM(i)	t-statistics	-5,8	-5,8	-5,3	-15,5	-13,8	-26,4
		t/t(between)	1,0	1,0	0,9	2,7	2,4	4,6
	BLK(i)	t-statistics	-3,5	-3,5	-3,3	-9,4	-8,3	-16,3
		t/t(between)	1,0	1,0	0,9	2,7	2,4	4,6
	ED(i)	t-statistics	9,3	9,3	1,0	24,8	21,9	42,9
		t/t(between)	1,0	1,0	0,1	2,7	2,4	4,6
			R-Between R-GLS	R-HT X1** RRBetween R-OLS	R-FEVD
	FEM(i)	t-statistics	-11,0	-6,6	-1,0	-28,1	-14,7	-49,9
		t/t(between)	1,9	1,2	0,2	4,9	2,5	8,7
	BLK(i)	t-statistics	-2,8	-3,6	-1,8	-7,4	-7,6	-12,8
		t/t(between)	0,8	1,0	0,5	2,1	2,2	3,7
	ED(i)	t-statistics	13,8	16,7	6,5	37,4	21,8	63,5
		t/t(between)	1,5	1,8	0,7	4,0	2,4	6,8

Table 4 :

 4 Contribution to R2 in Within regressions with forward selection.

	Variable	Rank	Partial R2	R2
	EXP	1	0.6505	0.6505
	EXP2	2	0.0060	0.6564
	SMSA	3	0.0005	0.6569
	UNION	4	0.0004	0.6574
	MS	5	0.0002	0.6576
	OCC	6	0.0002	0.6578
	WKS	7	0.0002	0.6580
	IND	8	0.0001	0.6581
	SOUTH	9	0.0000	0.6581

Table 5 :

 5 Contribution to R2 in Between regressions with forward selection.

	Variable	Rank	Partial R2	R2
	mED	1	0.2129	0.2129
	mFEM	2	0.1445	0.3575
	mSMSA	3	0.0581	0.4156
	mEXP	4	0.0420	0.4576
	mEXPSQ	5	0.0298	0.4874
	mBLK	6	0.0156	0.5030
	mOCC	7	0.0101	0.5131
	mUNION	8	0.0118	0.5249
	mWKS	9	0.0062	0.5311
	mIND	10	0.0051	0.5362
	mMS	11	0.0043	0.5405
	mSOUTH	12	0.0038	0.5443

Five of these thirteen estimates are also reported in[START_REF] Baltagi | A Companion to Econometric Analysis of Panel Data[END_REF], p.27-29: Within, Between, Restricted GLS, Restricted OLS and p.157, Restricted Hausman Taylor HT(X * * 1 ).