add eigvec field to the TwoVarCompVariateRotate type

src/VarianceComponentModels.jl

@@ -192,11 +192,12 @@ immutable TwoVarCompVariateRotate{T <: AbstractFloat, YT <: AbstractVecOrMat,
   Yrot::YT
   Xrot::XT
   eigval::Vector{T}
+  eigvec::Matrix{T}
   logdetV2::T
   # inner constructor
   function TwoVarCompVariateRotate(Yrot::AbstractVecOrMat{T},
-    Xrot::AbstractVecOrMat{T}, eigval::Vector{T}, logdetV2::T)
-    new(Yrot, Xrot, eigval, logdetV2)
+    Xrot::AbstractVecOrMat{T}, eigval::Vector{T}, eigvec::Matrix{T}, logdetV2::T)
+    new(Yrot, Xrot, eigval, eigvec, logdetV2)
   end
 end
 
@@ -209,10 +210,11 @@ function TwoVarCompVariateRotate(
   Yrot::AbstractVecOrMat,
   Xrot::AbstractVecOrMat,
   eigval::Vector,
+  eigvec::Matrix,
   logdetV2::Real)
 
   TwoVarCompVariateRotate{eltype(Yrot), typeof(Yrot), typeof(Xrot)}(Yrot, Xrot,
-    eigval, logdetV2)
+    eigval, eigvec, logdetV2)
 end
 
 """
@@ -225,6 +227,7 @@ function TwoVarCompVariateRotate{T <: AbstractFloat}(
   vcobs::VarianceComponentVariate{T, 2}
   )
 
+  zeroT = zero(T)
   # (generalized)-eigendecomposition of (V1, V2)
   if isa(vcobs.V[2], UniformScaling) ||
     (isdiag(vcobs.V[2]) && vecnorm(diag(vcobs.V[2]) - one(T)) < 1.0e-8)
@@ -239,12 +242,14 @@ function TwoVarCompVariateRotate{T <: AbstractFloat}(
     logdetV2 = convert(T, logdet(vcobs.V[2]))
   end
   # corect negative eigenvalues due to roundoff error
-  map!(x -> max(x, zero(T)), deval)::Vector{T}
+  @inbounds @simd for i in eachindex(deval)
+    deval[i] = deval[i] > zeroT ? deval[i] : zeroT
+  end
   # rotate responses
   Yrot = At_mul_B(U, vcobs.Y)
   Xrot = isempty(vcobs.X) ? Array{T}(size(Yrot, 1), 0) : At_mul_B(U, vcobs.X)
   # output
-  TwoVarCompVariateRotate(Yrot, Xrot, deval, logdetV2)
+  TwoVarCompVariateRotate(Yrot, Xrot, deval, U, logdetV2)
 end
 
 function TwoVarCompVariateRotate{T <: VarianceComponentVariate}(vcdata::Array{T})