diff --git a/DESCRIPTION b/DESCRIPTION
index ae16c55..ef629ca 100644
--- a/DESCRIPTION
+++ b/DESCRIPTION
@@ -1,6 +1,6 @@
 Package: rutilstimflutre
 Title: Timothee Flutre's personal R code
-Version: 0.152.0
+Version: 0.152.1
 Authors@R: c(
     person("Timothee", "Flutre", email="timothee.flutre@inra.fr", role=c("aut", "ctb", "cre")),
     person("Peter", "Carbonetto", role="aut", comment="function from varbvs example"),
diff --git a/R/GEMMA.R b/R/GEMMA.R
index e983e61..a784fa1 100644
--- a/R/GEMMA.R
+++ b/R/GEMMA.R
@@ -97,6 +97,8 @@ genoDoses2bimbam <- function(X=NULL, tX=NULL, alleles, file=NULL, format="mean")
 ##' ## test SNPs one by one with the univariate LMM
 ##' fit.u <- gemma(model="ulmm", y=modelA$Y[,1], X=X, snp.coords, alleles,
 ##'                W=modelA$W, out.dir=tempdir(), clean="all")
+##' fit.u$global.mean
+##' fit.u$pve
 ##' cor(modelA$a[modelA$gamma == 1], fit.u$tests$beta[modelA$gamma == 1])
 ##' cols <- rep("black",ncol(X)); cols[modelA$gamma==1] <- "red"
 ##' pvadj.AA <- qqplotPval(fit.u$tests$p_wald, col=cols, ctl.fdr.bh=TRUE,
@@ -302,6 +304,13 @@ gemma <- function(model="ulmm", y, X, snp.coords, alleles=NULL,
     output[["global.mean"]] <- c(beta.hat=as.numeric(tmp1[length(tmp1)]),
                                  se.beta.hat=as.numeric(tmp2[length(tmp2)]))
   }
+  idx <- grep("pve estimate in the null model", output[["log"]])
+  if(length(idx) == 1){
+    tmp1 <- strsplit(output$log[idx], " ")[[1]]
+    tmp2 <- strsplit(output$log[idx+1], " ")[[1]]
+    output[["pve"]] <- c(pve.hat=as.numeric(tmp1[length(tmp1)]),
+                         se.pve.hat=as.numeric(tmp2[length(tmp2)]))
+  }
   if(clean == "all")
     file.remove(f)
   if(model == "ulmm"){
diff --git a/man/gemma.Rd b/man/gemma.Rd
index 2d86aca..c722024 100644
--- a/man/gemma.Rd
+++ b/man/gemma.Rd
@@ -74,6 +74,8 @@ summary(abs(modelA$a[modelA$gamma == 1]))
 ## test SNPs one by one with the univariate LMM
 fit.u <- gemma(model="ulmm", y=modelA$Y[,1], X=X, snp.coords, alleles,
                W=modelA$W, out.dir=tempdir(), clean="all")
+fit.u$global.mean
+fit.u$pve
 cor(modelA$a[modelA$gamma == 1], fit.u$tests$beta[modelA$gamma == 1])
 cols <- rep("black",ncol(X)); cols[modelA$gamma==1] <- "red"
 pvadj.AA <- qqplotPval(fit.u$tests$p_wald, col=cols, ctl.fdr.bh=TRUE,