Calculate repeatability

calculate_repeatability.r

@@ -0,0 +1,102 @@
+# calculate_repeatability.r
+
+library(lme4)
+library(heritability)
+
+
+
+# Input long format data
+path.rgb.long <- ("./df_rgb_long.csv")
+path.hsi.long <- ("./df_hsi_long.csv")
+path.rpt      <- ("./rpt.csv")
+
+factor.col <- c('EXP.ID', 
+                'POT_BARCODE',
+                'TREATMENT', 
+                'VARIETY', 
+                'GROWTH_STAGE',
+                'View',
+                'frame_nr')
+
+
+# Load RGB and HSI data in long format
+message("Loading data...")
+df.rgb <- read.csv(path.rgb.long)
+df.hsi <- read.csv(path.hsi.long)
+
+# Combine RGB and HSI data
+df <- rbind(df.rgb, df.hsi)
+
+# Set factor columns
+df[, factor.col] <- lapply(df[, factor.col], factor)
+
+# Create a dataframe to record repeatability
+message("Finding unique conditions to evaluate repeatability...")
+df.rpt <- unique(df[, c('TREATMENT', 
+                        'GROWTH_STAGE', 
+                        'View', 
+                        'frame_nr', 
+                        'variable')])
+df.rpt <- df.rpt[order(df.rpt$TREATMENT, 
+                       df.rpt$GROWTH_STAGE, 
+                       df.rpt$View, 
+                       df.rpt$frame_nr, 
+                       df.rpt$variable), ]
+df.rpt$repeatability <- NA
+df.rpt$gen.variance  <- NA
+df.rpt$res.variance  <- NA
+
+# Calculate repeatability
+message('Calculating repeatability...')
+for(i in 1:nrow(df.rpt)){
+    row <- df.rpt[i,]
+    # Subset data by condition by which repeatability is measured
+    message(paste(row$TREATMENT, 
+                  row$GROWTH_STAGE, 
+                  row$View, 
+                  row$frame_nr, 
+                  row$variable, 
+                  " (", i, "/", nrow(df.rpt), ")" sep=" "))
+    df.temp <- df[which(df$TREATMENT==row$TREATMENT & 
+                        df$GROWTH_STAGE==row$GROWTH_STAGE & 
+                        df$View==row$View & 
+                        df$frame_nr==row$frame_nr & 
+                        df$variable==row$variable), ]
+
+    # Count NA values and continue if excessive
+    num.na <- sum(is.na(df.temp$value))
+    if (num.na/nrow(df.temp) > 0.3) next
+
+    # Select valid data
+    df.temp <- df.temp[!is.na(df.temp$value) & is.finite(df.temp$value), ]
+
+    # Evaluate repeatability (line repeatability=True)
+    rpt <- repeatability(df.temp$value, df.temp$VARIETY, line.repeatability=T,
+                         covariates.frame = data.frame())
+
+    # Save data
+    df.rpt$repeatability[i] <- rpt$repeatability
+    df.rpt$gen.variance[i]  <- rpt$gen.variance
+    df.rpt$res.variance[i]  <- rpt$res.variance
+}
+
+
+
+
+# Export repeatability table
+message('\nExporting repeatability measurements...')
+write.csv(df.rpt, path.rpt, row.names=F)
+
+
+
+
+# factor(EXP.ID, frame_nr, VARIETY, TREATMENT)
+# value ~ 1|EXP.ID + VARIETY + TREATMENT + GROWTH_STAGE
+# value ~ 1|EXP.ID + VARIETY + TREATMENT                Seth
+# value ~ 1|EXP.ID + VARIETY                            Current
+
+
+
+
+
+# EOF