diff --git a/NAMESPACE b/NAMESPACE
index 51527ed..2e4c98a 100644
--- a/NAMESPACE
+++ b/NAMESPACE
@@ -3,7 +3,12 @@
 export(comb.indices)
 export(gen.advance)
 export(gen.varcov)
+export(meanPerformance)
 export(phen.varcov)
 export(rcomb.indices)
 export(weight.mat)
+importFrom(stats,aggregate)
+importFrom(stats,anova)
+importFrom(stats,lm)
+importFrom(stats,qt)
 importFrom(utils,combn)
diff --git a/R/meanPerformance.R b/R/meanPerformance.R
new file mode 100644
index 0000000..09fdf71
--- /dev/null
+++ b/R/meanPerformance.R
@@ -0,0 +1,95 @@
+#' @title Mean performance of phenotypic data
+#'
+#' @param data data for analysis
+#' @param genotypes genotypes vector
+#' @param replications replication vector
+#'
+#' @return Dataframe of mean performance analysis
+#' @export
+#' @importFrom stats aggregate anova lm qt
+#' @examples
+#' meanPerformance(data = seldata[, 3:9], genotypes = seldata[, 2], replications = seldata[, 1])
+#'
+#'
+
+meanPerformance <- function(data, genotypes, replications){
+  # Calculate average results from replicaions
+  meanData <- function(data, genotypes){
+    odr <- unique(genotypes)
+    genotype <- list(factor(genotypes, ordered(odr)))
+    mean <- aggregate(data, genotype, mean)
+    matrix <- as.matrix(mean[,-1])
+    matrix <- round(matrix, 4)
+    result <- data.frame("Genotypes" = odr, matrix)
+    return(data.frame(result))
+  }
+  # Calculate the performance
+  performance <- function (data, treat, rep)
+  {
+    convert <- function(data) {
+      data <- as.data.frame(sapply(data, as.numeric))
+      data <- as.list(data)
+      return(data)
+    }
+    analysis <- function(data, treat, rep) {
+      data <- as.numeric(data)
+      genotype <- as.factor(treat)
+      replication <- as.factor(rep)
+      r <- nlevels(replication)
+      model <- lm(data ~ replication + genotype)
+      anova.model <- anova(model)
+      MinMax <- meanData(data = data, genotypes = treat)
+      Maxi <- round(max(MinMax[, -1]), 4)
+      Mini <- round(min(MinMax[, -1]), 4)
+      rang <- paste(round(Maxi, 4), " - ", round(Mini, 4))
+      GM <- round(mean(data), 4)
+      EMS <- round(anova.model[3, 3],4)
+      SD <- round(sqrt(EMS),4)
+      SEm <- round(sqrt(EMS/r), 4)
+      CV<- round(SD/GM*100, 4)
+      CD5 <- round(sqrt(EMS/r) * sqrt(2) * abs(qt(0.025, anova.model[3, 1])), 4)
+      if (anova.model[2, 5] > 0.05) {
+        CD5 <- paste(round(sqrt(EMS/r) * sqrt(2) * abs(qt(0.025,
+                                                          anova.model[3, 1])), 4), "NS")
+      }
+      CD1 <- round(sqrt(EMS/r) * sqrt(2) * abs(qt(0.005, anova.model[3,
+                                                                     1])), 4)
+      if (anova.model[2, 5] > 0.01) {
+        CD1 <- paste(round(sqrt(EMS/r) * sqrt(2) * abs(qt(0.005,
+                                                          anova.model[3, 1])), 4), "NS")
+      }
+      GV <- round((anova.model[2, 3] - EMS)/r, 4)
+      PV <- round(GV + EMS, 4)
+      hs <- round((GV/PV), 4)
+      matri <- matrix(data = c(Mini, Maxi, GM, CV, SEm, CD5, CD1, hs, hs*100),
+                      dimnames = list(c("Min", "Max", "Grand Mean", "CV %",
+                                        "SEm",
+                                        "CD 5%",
+                                        "CD 1%",
+                                        "Heritability",
+                                        "Heritability in %")),
+                      nrow = 9)
+      table1 <- as.data.frame(matri, useNa = F)
+      my.list <- list(table1)
+      return(my.list)
+    }
+    fiftn <- convert(data)
+    colnumber <- ncol(data)
+    output <- list()
+    for (j in 1:colnumber) {
+      output[[j]] <- analysis(fiftn[[j]], treat,
+                              rep)
+    }
+    names(output) <- names(data)
+    return(output)
+  }
+  meandf <- meanData(data = data, genotypes = genotypes)
+  performanceList <- performance(data = data, treat = genotypes, rep = replications)
+  performanceBind <- cbind.data.frame(performanceList, row.names = NULL)
+  Genotypes <- matrix(c("Min", "Max", "GM", "CV (%)", "SEm", "CD 5%", "CD 1%", "Heritability", "Heritability(%)"),
+                      ncol = 1)
+  FinalPerformanceBind <- cbind.data.frame(Genotypes, performanceBind, row.names = NULL)
+  colnames(FinalPerformanceBind) <- colnames(meandf)
+  result <- rbind.data.frame(meandf, FinalPerformanceBind)
+  return(result)
+}
diff --git a/man/meanPerformance.Rd b/man/meanPerformance.Rd
new file mode 100644
index 0000000..14a5051
--- /dev/null
+++ b/man/meanPerformance.Rd
@@ -0,0 +1,26 @@
+% Generated by roxygen2: do not edit by hand
+% Please edit documentation in R/meanPerformance.R
+\name{meanPerformance}
+\alias{meanPerformance}
+\title{Mean performance of phenotypic data}
+\usage{
+meanPerformance(data, genotypes, replications)
+}
+\arguments{
+\item{data}{data for analysis}
+
+\item{genotypes}{genotypes vector}
+
+\item{replications}{replication vector}
+}
+\value{
+Dataframe of mean performance analysis
+}
+\description{
+Mean performance of phenotypic data
+}
+\examples{
+meanPerformance(data = seldata[, 3:9], genotypes = seldata[, 2], replications = seldata[, 1])
+
+
+}