PROBLEM STATEMENT
I have a code to treat agronomic data. I have a data.frame called yield
with different information (dry.weight / infected pod / healthy pod). These data are collected in the field every month. I would like to have the total yield per year and per tree (plotnb).
Another important thing I'm seeking is the total number of pods produced (either diseased or healthy) the problem I have is that some pods are preleved to make some disease resistance test and therefore are not counted as either healthy or diseased.
It is a routine code and (with whole data) it takes around 1.5 minutes to calculate and I find it too much for what I'm calculating so I'm sure there is a cleverest way to do it.
I'm seeking a code review here to see what you guys think of my code and if there is a way to optimize it and I will take any advice you could give me to improve it!
REPRODUCIBLE DATA
yield <- structure(list(plotnb = c(49L, 49L, 49L, 49L, 49L, 49L, 49L,
49L, 49L, 49L, 49L, 89L, 158L, 158L, 158L, 158L, 158L, 158L,
158L, 158L, 158L, 159L, 159L, 249L, 249L, 249L, 318L, 318L, 318L,
326L, 326L, 326L, 326L, 326L, 349L, 349L, 408L, 421L, 421L, 421L,
421L, 421L, 423L, 423L, 423L, 424L, 424L, 424L, 424L, 424L, 424L,
506L, 506L, 506L, 562L, 562L, 562L, 562L, 562L, 562L, 562L, 562L,
562L, 562L, 562L, 562L, 562L, 649L, 649L, 747L, 747L, 747L, 747L,
747L, 747L, 798L, 866L, 866L, 866L, 866L, 866L, 930L, 930L, 930L,
930L, 930L, 930L, 930L, 930L, 930L, 930L, 930L, 963L, 963L, 963L,
963L, 963L, 963L, 963L, 963L, 1016L, 1016L, 1016L, 1016L, 1016L,
1016L, 1016L, 1016L, 1066L, 1066L, 1102L, 1102L, 1102L, 1102L,
1102L, 1185L, 1185L, 1185L, 1185L, 1185L, 1185L, 1185L, 1185L,
1185L, 1185L, 1185L, 1185L, 1186L, 1186L, 1186L, 1186L, 1186L,
1186L, 1186L, 1194L, 1194L, 1194L, 1194L, 1435L, 1531L, 1531L,
1531L, 1531L, 1531L, 1531L, 1531L, 1547L, 1559L, 1559L, 1559L,
1559L, 1559L),
dry.weight = c(24L, 116L, 52L, 30L, 142L, 40L,
34L, 10L, 52L, 26L, 44L, 48L, 10L, 56L, 40L, 38L, 46L, 36L, 14L,
24L, 130L, 34L, 24L, 56L, 30L, 28L, 52L, 386L, 46L, 46L, 16L,
28L, 32L, 28L, 22L, 28L, 22L, 58L, 14L, 40L, 14L, 96L, 142L,
114L, 46L, 34L, 46L, 114L, 130L, 38L, 134L, 44L, 42L, 26L, 34L,
42L, 18L, 10L, 40L, 102L, 56L, 24L, 12L, 44L, 46L, 18L, 30L,
52L, 58L, 52L, 4L, 64L, 14L, 74L, 206L, 30L, 108L, 20L, 46L,
6L, 40L, 46L, 28L, 32L, 102L, 68L, 58L, 48L, 32L, 74L, 32L, 114L,
58L, 32L, 28L, 48L, 6L, 32L, 26L, 64L, 108L, 34L, 46L, 84L, 28L,
84L, 34L, 88L, 20L, 46L, 66L, 152L, 164L, 48L, 84L, 470L, 70L,
42L, 294L, 110L, 174L, 126L, 54L, 872L, 48L, 312L, 62L, 162L,
44L, 46L, 90L, 34L, 228L, 188L, 78L, 406L, 170L, 168L, 36L, 36L,
76L, 24L, 30L, 58L, 82L, 124L, 32L, 76L, 36L, 88L, 94L, 26L),
healthy.pod = c(1L, 2L, 1L, 1L, 5L, 2L, 1L, 1L, 1L, 1L, 2L,
1L, 2L, 3L, 2L, 1L, 1L, 1L, 2L, 1L, 6L, 1L, 1L, 1L, 1L, 1L,
2L, 17L, 2L, 5L, 2L, 1L, 1L, 4L, 1L, 1L, 1L, 1L, 1L, 1L,
1L, 5L, 14L, 8L, 2L, 1L, 2L, 6L, 5L, 2L, 8L, 2L, 1L, 1L,
2L, 1L, 1L, 1L, 1L, 4L, 1L, 1L, 1L, 2L, 1L, 1L, 1L, 1L, 2L,
1L, 1L, 1L, 1L, 2L, 7L, 1L, 3L, 1L, 2L, 1L, 1L, 1L, 1L, 1L,
2L, 2L, 1L, 1L, 1L, 2L, 1L, 3L, 4L, 1L, 2L, 2L, 1L, 3L, 1L,
2L, 8L, 1L, 2L, 4L, 1L, 2L, 2L, 3L, 1L, 1L, 1L, 3L, 5L, 2L,
3L, 12L, 2L, 2L, 7L, 2L, 4L, 3L, 1L, 19L, 1L, 6L, 1L, 4L,
1L, 1L, 2L, 1L, 5L, 3L, 2L, 13L, 3L, 4L, 1L, 1L, 2L, 2L,
1L, 4L, 2L, 3L, 1L, 1L, 1L, 3L, 2L, 2L),
infected.pods = c(0L,
0L, 0L, 4L, 2L, 2L, 0L, 0L, 0L, 0L, 2L, 0L, NA, 1L, 0L, 0L,
0L, 0L, 0L, NA, 0L, 0L, NA, 0L, 0L, 0L, 2L, 0L, 0L, 1L, 0L,
0L, 0L, 3L, 0L, 2L, 0L, 2L, NA, NA, 0L, 0L, 2L, 0L, 0L, 0L,
2L, 1L, 2L, 3L, 4L, 0L, 2L, NA, 0L, 0L, 0L, 0L, 0L, 0L, 0L,
0L, 0L, NA, 0L, 0L, NA, 0L, 0L, NA, NA, 0L, 0L, 0L, 1L, 0L,
0L, 0L, 0L, 0L, NA, NA, 2L, 2L, 0L, 1L, 0L, 0L, 0L, 1L, 0L,
0L, NA, 2L, 0L, 2L, 0L, 1L, 0L, NA, 0L, 0L, 0L, 1L, NA, NA,
0L, 0L, 0L, 0L, 0L, NA, 0L, 1L, 2L, 6L, 0L, 0L, 0L, 0L, NA,
1L, 0L, 0L, 0L, 0L, 1L, 1L, NA, 0L, 0L, 0L, 0L, 0L, 0L, NA,
0L, 0L, 0L, 0L, 0L, 0L, NA, 0L, 0L, 0L, NA, 0L, 0L, 0L, 2L,
NA),
date = structure(c(29L, 35L, 37L, 5L, 25L, 9L, 16L,
13L, 33L, 7L, 11L, 11L, 8L, 3L, 19L, 7L, 12L, 17L, 4L, 36L,
16L, 33L, 36L, 23L, 35L, 12L, 5L, 16L, 13L, 33L, 31L, 24L,
9L, 37L, 17L, 9L, 35L, 11L, 36L, 8L, 33L, 7L, 25L, 29L, 28L,
28L, 16L, 36L, 5L, 13L, 9L, 27L, 9L, 8L, 25L, 1L, 28L, 11L,
35L, 16L, 22L, 5L, 29L, 36L, 34L, 31L, 8L, 5L, 17L, 36L,
38L, 16L, 3L, 13L, 9L, 24L, 37L, 9L, 28L, 35L, 36L, 8L, 5L,
23L, 28L, 9L, 34L, 37L, 33L, 12L, 15L, 35L, 36L, 15L, 17L,
16L, 9L, 8L, 5L, 38L, 27L, 13L, 35L, 23L, 8L, 36L, 17L, 12L,
17L, 24L, 24L, 36L, 35L, 5L, 16L, 16L, 33L, 25L, 29L, 28L,
8L, 17L, 13L, 24L, 5L, 27L, 9L, 25L, 8L, 37L, 13L, 35L, 29L,
28L, 25L, 36L, 28L, 29L, 30L, 30L, 17L, 26L, 39L, 37L, 24L,
32L, 39L, 20L, 26L, 30L, 21L, 39L),
.Label = c("02/09/2015", "03/08/2015", "04/07/2016", "04/08/2016", "04/08/2017", "04/09/2016",
"05/05/2016", "05/10/2017", "06/07/2017", "06/10/2017", "07/04/2016",
"07/04/2017", "07/06/2017", "07/07/2015", "07/09/2016", "07/09/2017",
"07/10/2016", "08/01/2018", "08/06/2016", "08/06/2017", "08/08/2017",
"08/10/2015", "09/05/2017", "09/12/2016", "10/03/2016", "10/05/2017",
"10/11/2016", "11/01/2016", "11/02/2016", "11/09/2017", "11/11/2015",
"11/11/2016", "12/01/2017", "14/12/2015", "16/03/2017", "16/11/2017",
"17/02/2017", "18/12/2017", "20/11/2017"), class = "factor")), .Names = c("plotnb",
"dry.weight", "healthy.pod", "infected.pods", "date"), row.names = c(286L,
287L, 288L, 289L, 290L, 291L, 292L, 293L, 294L, 295L, 296L, 503L,
924L, 925L, 926L, 927L, 928L, 929L, 930L, 931L, 932L, 933L, 934L,
1365L, 1366L, 1367L, 1790L, 1791L, 1792L, 1846L, 1847L, 1848L,
1849L, 1850L, 1981L, 1982L, 2366L, 2450L, 2451L, 2452L, 2453L,
2454L, 2458L, 2459L, 2460L, 2461L, 2462L, 2463L, 2464L, 2465L,
2466L, 2962L, 2963L, 2964L, 3212L, 3213L, 3214L, 3215L, 3216L,
3217L, 3218L, 3219L, 3220L, 3221L, 3222L, 3223L, 3224L, 3531L,
3532L, 3971L, 3972L, 3973L, 3974L, 3975L, 3976L, 4166L, 4387L,
4388L, 4389L, 4390L, 4391L, 4605L, 4606L, 4607L, 4608L, 4609L,
4610L, 4611L, 4612L, 4613L, 4614L, 4615L, 4747L, 4748L, 4749L,
4750L, 4751L, 4752L, 4753L, 4754L, 5030L, 5031L, 5032L, 5033L,
5034L, 5035L, 5036L, 5037L, 5252L, 5253L, 5411L, 5412L, 5413L,
5414L, 5415L, 5761L, 5762L, 5763L, 5764L, 5765L, 5766L, 5767L,
5768L, 5769L, 5770L, 5771L, 5772L, 5773L, 5774L, 5775L, 5776L,
5777L, 5778L, 5779L, 5794L, 5795L, 5796L, 5797L, 6620L, 6807L,
6808L, 6809L, 6810L, 6811L, 6812L, 6813L, 6840L, 6854L, 6855L,
6856L, 6857L, 6858L), class = "data.frame")
monilia <- structure(list(plotnb = structure(c(24L, 24L, 24L, 156L, 162L,
162L, 179L, 218L, 219L, 219L, 237L, 237L, 237L, 332L, 332L, 385L,
385L, 385L), .Label = c("1", "10", "100", "101", "102", "103",
"106", "107", "1073", "1074", "1078", "1079", "108", "1082",
"1086", "1088", "109", "1091", "1097", "1098", "1099", "11",
"1101", "1102", "1104", "1105", "1106", "1107", "1108", "1109",
"111", "1112", "1116", "1118", "112", "1124", "1127", "1128",
"1129", "1133", "1134", "1136", "1137", "1139", "1142", "1145",
"1146", "1148", "115", "1151", "1152", "1153", "1155", "1157",
"116", "1173", "118", "1180", "119", "1201", "121", "1242", "1243",
"1248", "1259", "126", "1260", "1280", "1281", "1290", "1299",
"13", "1302", "1318", "1334", "1347", "1365", "1375", "14", "1403",
"1407", "1408", "141", "1412", "1445", "1446", "1447", "1451",
"1452", "1453", "1455", "1457", "1467", "1472", "1476", "1483",
"1492", "1519", "1524", "1525", "160", "172", "179", "18", "182",
"183", "192", "2", "203", "21", "22", "220", "23", "235", "240",
"247", "257", "258", "259", "26", "260", "261", "262", "264",
"27", "271", "273", "274", "275", "276", "277", "278", "279",
"28", "281", "283", "288", "290", "291", "292", "293", "295",
"297", "298", "301", "302", "303", "305", "306", "307", "309",
"31", "310", "313", "314", "318", "320", "321", "323", "324",
"325", "326", "33", "331", "332", "334", "335", "336", "337",
"34", "340", "341", "342", "343", "344", "346", "347", "348",
"349", "350", "352", "355", "356", "357", "359", "360", "363",
"364", "367", "368", "369", "372", "373", "375", "376", "377",
"379", "38", "380", "381", "383", "386", "389", "391", "394",
"398", "399", "4", "40", "400", "406", "41", "410", "411", "413",
"414", "419", "423", "424", "430", "431", "433", "436", "437",
"44", "441", "442", "443", "447", "453", "457", "46", "461",
"470", "479", "486", "49", "490", "491", "497", "5", "50", "503",
"52", "521", "524", "530", "533", "536", "539", "542", "547",
"551", "552", "554", "556", "558", "561", "564", "568", "57",
"574", "577", "579", "58", "580", "581", "582", "587", "588",
"589", "590", "593", "597", "598", "599", "602", "604", "606",
"607", "609", "61", "611", "613", "617", "62", "63", "637", "65",
"655", "657", "66", "661", "662", "664", "666", "668", "67",
"671", "683", "684", "685", "686", "688", "696", "698", "7",
"702", "703", "704", "706", "71", "710", "711", "712", "717",
"718", "721", "722", "724", "726", "733", "734", "735", "737",
"739", "74", "740", "742", "743", "746", "747", "748", "755",
"758", "76", "761", "762", "77", "773", "774", "777", "778",
"78", "781", "783", "786", "789", "796", "797", "80", "801",
"803", "804", "807", "808", "815", "816", "817", "819", "82",
"820", "822", "823", "824", "826", "827", "828", "83", "830",
"831", "836", "837", "838", "84", "844", "852", "853", "856",
"858", "86", "863", "864", "865", "866", "87", "872", "876",
"877", "880", "92", "94", "98", "99", "RDT53"), class = "factor"),
Fecha.Calificacion = structure(c(3L, 3L, 3L, 6L, 5L, 10L,
5L, 9L, 8L, 8L, 14L, 14L, 14L, 16L, 16L, 12L, 12L, 12L), .Label = c("04/07/2017",
"04/10/2017", "12/09/2017", "13/11/2017", "15/08/2017", "16/05/2017",
"17/05/2017", "18/05/2017", "20/06/2017", "20/09/2017", "20/12/2017",
"23/08/2017", "23/10/2017", "25/07/2017", "26/06/2017", "28/06/2017",
"28/11/2017", "29/11/2017", "30/05/2017", "31/05/2017"), class = "factor")), .Names = c("plotnb",
"Fecha.Calificacion"), row.names = c(59L, 60L, 61L, 400L, 412L,
413L, 456L, 552L, 553L, 554L, 591L, 592L, 593L, 768L, 769L, 907L,
908L, 909L), class = "data.frame")
phytophtora <-structure(list(plotnb = structure(c(17L, 17L, 17L, 17L, 80L,
80L, 80L), .Label = c("1072", "1073", "1074", "1075", "1078",
"1082", "1086", "1087", "1088", "1091", "1093", "1097", "1098",
"1099", "1100", "1101", "1102", "1104", "1106", "1108", "1109",
"1112", "1116", "1122", "1127", "1128", "1129", "1130", "1131",
"1136", "1138", "1139", "1141", "1142", "1143", "1144", "1146",
"1148", "1150", "1151", "1153", "1154", "1157", "1159", "375",
"777", "778", "779", "781", "783", "788", "796", "799", "801",
"803", "804", "807", "809", "812", "816", "819", "820", "823",
"824", "827", "828", "836", "837", "838", "842", "843", "845",
"846", "856", "858", "859", "861", "863", "864", "866", "867",
"869", "871", "872", "875", "877", "RDT10"), class = "factor"),
Fecha.Calificacion = structure(c(3L, 7L, 1L, 2L, 7L, 7L,
7L), .Label = c("08/05/2017", "10/04/2017", "14/08/2017",
"15/09/2017", "16/11/2017", "25/01/2018", "29/06/2017"), class = "factor")), .Names = c("plotnb",
"Fecha.Calificacion"), row.names = c(36L, 37L, 38L, 39L, 170L,
171L, 172L), class = "data.frame")
CODE
## using a sys.time() to check time taken.
start.time <- Sys.time()
## creating a column year to better group data.
yield$year <- "2015-2016"
yield$year[as.POSIXct(as.character(yield$date),format="%d/%m/%Y") >= "2017-07-01"] <- "2017-2018"
yield$year[as.POSIXct(as.character(yield$date),format="%d/%m/%Y") >= "2016-07-01" &
as.POSIXct(as.character(yield$date),format="%d/%m/%Y") < "2017-07-01"] <- "2016-2017"
yield$year <- as.factor(yield$year)
## grouping data and calculating Total pod/ Total Healthy for each year.
yield.year <- c()
for (parcelle in unique(yield$plotnb)){
for (year in levels(yield$year)){
yield.subset <- yield[yield$plotnb==parcelle & yield$year==year,]
TotalSane= sum(yield.subset$healthy.pod)
TotalInfected= sum(yield.subset$infected.pods)
TotalPreleved=0
if (year=="2016-2017"){
monilia.subset= monilia[monilia$plotnb==parcelle & as.POSIXct(as.character(monilia$Fecha.inoculacion),format="%d/%m/%Y") < "2017-07-01" ,]
phytophtora.subset= phytophtora[phytophtora$plotnb==parcelle & as.POSIXct(as.character(phytophtora$Fecha.inoculacion),format="%d/%m/%Y") < "2017-07-01" ,]
TotalPreleved = nrow(monilia.subset)+ nrow(phytophtora.subset)
} else if (year=="2017-2018"){
monilia.subset= monilia[monilia$plotnb==parcelle & as.POSIXct(as.character(monilia$Fecha.inoculacion),format="%d/%m/%Y") >= "2017-07-01" ,]
phytophtora.subset= phytophtora[phytophtora$plotnb==parcelle & as.POSIXct(as.character(phytophtora$Fecha.inoculacion),format="%d/%m/%Y") >= "2017-07-01" ,]
TotalPreleved = nrow(monilia.subset)+ nrow(phytophtora.subset)
}
TotalPod= TotalSane + TotalInfected + TotalPreleved
TotalWeight= sum(yield.subset$dry.weight)
yield.year <- rbind(yield.year,c(parcelle,TotalSane,TotalInfected,TotalPreleved,TotalPod,TotalWeight,year))
}
}
## formating
yield.year <- as.data.frame(yield.year)
colnames(yield.year) <- c("plotnb","TotalSane","TotalInfected","TotalPreleved","TotalPod","TotalWeight","year")
## calculating other data
yield.year[,"potential_yield"] <- as.numeric(as.vector(yield.year$TotalPod))*as.numeric(as.vector(yield.year$TotalWeight))/as.numeric(as.vector(yield.year$TotalSane))
yield.year[,"weight_per_pod"] <- as.numeric(as.vector(yield.year$TotalWeight))/as.numeric(as.vector(yield.year$TotalSane))
yield.year[,"TotalWeight"] <- as.numeric(as.character(yield.year$TotalWeight))
yield.year[,"TotalPod"] <- as.numeric(as.character(yield.year$TotalPod))
## cleaning
rm(list=c("monilia.subset","phytophtora.subset","yield.subset",
"parcelle","TotalInfected","TotalPod","TotalPreleved","TotalSane","TotalWeight","year"))
## getting time.taken
end.time <- Sys.time()
time.taken <- end.time - start.time
time.taken
TotalPreleved = nrow(monilia.subset)+ nrow(phytophtora.subset)
as the second data set is specified later if at all. So I skipped this task completely ou in my answer. Please revise your question how to summarisemonilia
andphytophtora
. Concrete how to sum up forplotnb == 49
theFecha.Calificacion
. By the way it is three times the same date. \$\endgroup\$plotnb
andYear
and add this information. \$\endgroup\$monilia.subset
is not3
foryear=="2017-2018"
andplotnb==49
? \$\endgroup\$