Расчет энергии впуска с помощью петель - PullRequest
Новая
       4

Расчет энергии впуска с помощью петель

0 голосов
/ 23 апреля 2020

Я пытаюсь запустить старый сценарий коллег и надеюсь, что кто-нибудь может помочь мне сообщить, что именно он делал во время этого сегмента кода. Ранее в сценарии мы рассчитывали норму потребления для нескольких видов добычи, и теперь кажется, что мы группируем их по уникальным местам. Раздел кода после этого требует, чтобы было 41 строка (1 строка для каждого уникального местоположения в полном наборе данных). Я считаю, что код размещает данные на основе широты, а затем добавляет столбец «альфа». Основной вопрос, который у меня возник, заключается в том, что вычисляет эта строка: x= x + d$Intakerate_kjday[j]*d$alpha[j]. Для локаций, у которых было несколько предметов добычи (profit.fall.38.4.959), этот код суммирует «capturerate_kjday» и «alpha», а затем умножает их вместе? Когда код выполняется, я получаю ошибку Error в 

 `[<-.data.frame`(`*tmp*`, k, , value = c("2", "Bishop's Head",  : replacement has 6 items, need 7

. Я был бы очень признателен за понимание того, что он пытался вычислить, и возможного обходного пути. Спасибо. 

dput(profit)
structure(list(Sample.ID = structure(c(5L, 19L, 27L, 28L, 30L, 
38L, 12L, 62L, 49L, 29L, 25L, 17L, 61L, 67L, 27L, 26L, 32L, 9L, 
47L, 45L, 5L, 26L, 27L, 44L, 45L, 4L, 1L, 43L, 19L, 35L), .Label = c("Barren Island Mud 1", 
"BH High 1", "BH High 2", "BH Low 1", "BH Low 2", "BH Low 3", 
"BHH 1 C", "BHH 2 E", "BHL 1 E", "BHL 2", "BHL 3 (B)", "BHM 1 C", 
"BI High 1", "BI Low 1", "BI Low 2C", "BI Low 3", "BI Mud", "BIHI High B", 
"BIL1 (low) E", "BIL1E", "BIL2 E", "BIL2E", "BW Fresh 1", "BW Fresh 2", 
"BW High 1", "BW High 2", "BW High 5", "BW Low 3", "BW Money Stump", 
"BW Mud 1", "BW SAV 1", "BW SAV 2", "BWH 1 D", "BWH 2", "BWH 3", 
"BWH 5", "BWL 1", "BWL 2", "BWL 3", "BWM 1", "BWMS D", "EN High 2", 
"EN High 4", "EN High 5", "EN Low 1", "EN Low 2", "EN Mud 2", 
"ENH3 A High", "ENH4 A High", "ENH5 A High", "ENL1 Low E", "ENM1 A Mud", 
"ENS1 SAV", "ENS2 SAV 2C", "ENS3 SAV 3E", "High 3C", "MWP 29 Low 1", 
"MWP 30 Mud 1", "MWP 31 Low 2", "MWP 32 Mud 2", "MWP 33 Low 3", 
"MWP 34 Low 4", "PWRC Fresh", "WP 27 HM-MARC", "WP 28 HM-MARC", 
"WP 30 IT MARE", "WP29 LM-MARC"), class = "factor"), Season = structure(c(2L, 
3L, 2L, 2L, 2L, 3L, 3L, 2L, 3L, 2L, 2L, 2L, 2L, 3L, 2L, 2L, 2L, 
3L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 2L, 4L, 2L, 3L, 3L), .Label = c("", 
"Fall", "Spring", "Spring?"), class = "factor"), Refuge = structure(c(3L, 
2L, 5L, 5L, 5L, 5L, 4L, 7L, 6L, 5L, 5L, 2L, 7L, 7L, 5L, 5L, 5L, 
4L, 6L, 6L, 3L, 5L, 5L, 6L, 6L, 3L, 2L, 6L, 2L, 5L), .Label = c("", 
"Barren Island", "Bishop's Head", "Bishops Head", "Blackwater", 
"Eastern Neck", "Martin", "PWRC"), class = "factor"), Habitat.Type = structure(c(3L, 
3L, 2L, 3L, 4L, 3L, 4L, 3L, 2L, 3L, 2L, 4L, 3L, 3L, 2L, 2L, 5L, 
3L, 4L, 3L, 3L, 2L, 2L, 2L, 3L, 3L, 4L, 2L, 3L, 2L), .Label = c("Fresh", 
"High", "Low", "Mud", "SAV"), class = "factor"), Longitude = c(-76.03896, 
-76.26205, -76.05714, -76.06332, -76.14641, -76.23522, -76.03869, 
-75.99733, -76.21661, -76.23491, -76.22003, -76.26163, -75.99354, 
-76.01407, -76.05714, -76.01762, -76.10363, -76.04883, -76.21547, 
-76.23986, -76.03896, -76.01762, -76.05714, -76.2181, -76.23986, 
-76.04883, -76.26163, -76.21661, -76.26205, -76.0235), Latitude = c(38.22447, 
38.33905, 38.40959, 38.39708, 38.41795, 38.43055, 38.23255, 37.99369, 
39.03264, 38.43141, 38.41026, 38.33606, 37.98833, 38.01108, 38.40959, 
38.41913, 38.40351, 38.22694, 39.04036, 39.02677, 38.22447, 38.41913, 
38.40959, 39.03887, 39.02677, 38.22694, 38.33606, 39.03264, 38.33905, 
38.39138), Prey = structure(c(11L, 41L, 35L, 30L, 41L, 41L, 41L, 
3L, 18L, 31L, 40L, 9L, 41L, 38L, 30L, 13L, 35L, 41L, 20L, 27L, 
4L, 40L, 13L, 35L, 41L, 5L, 5L, 15L, 22L, 20L), .Label = c("Hydrobia", 
"Hydrobia genus", "Hydrobia sp.", "Hydrobia spp", "Melampus bidentatus", 
"Ruppia (maritima or rostellata)", "Ruppia genus", "Ruppia maritima", 
"Schoenoplectus pungens", "Schoenoplectus robustus", "Schoenoplectus spp", 
"Schoenoplectus spp.", "Scirpus acutus", "Scirpus acutus?", "Scirpus americanus", 
"Scirpus fluviatilis", "Scirpus genus", "Scirpus genus 1", "Scirpus genus 1?", 
"Scirpus genus 2", "Scirpus genus 3", "Scirpus genus?", "Scirpus heterochaetus", 
"Scirpus meterochaetus", "Scirpus mevadensis", "Scirpus olney?", 
"Scirpus olneyi", "Scirpus paludosis", "Scirpus paludosus", "Scirpus robustus", 
"Scirpus robustus?", "Scirpus species", "Scirpus subterminalis", 
"Scirpus subtermiralis", "Scirpus validus", "Spartina alterniflora", 
"Spartina genus", "Spartina genus?", "Spartina patens", "Spartina pectinata", 
"Zannichallia palustris"), class = "factor"), Density = c(2.36e-05, 
0.000101477, 0.000335244, 1.17e-05, 1.91e-06, 2.8e-06, 1.72e-05, 
1.34e-05, 2.71e-05, 0.000107843, 2.16e-06, 4.46e-06, 1.22e-05, 
6.61e-05, 0.000263052, 3.91e-05, 0.00034925, 3.69e-06, 8.02e-06, 
2.04e-05, 2.9e-05, 2.05e-05, 0.000564046, 0.001912535, 2.04e-05, 
0.001117905, 0.00255132, 9.03e-05, 4.23e-05, 0.000248282), Intakerate_kcals = c(-3.5399430250046e-07, 
7.6382794280604e-14, -5.02872205332896e-06, -1.7549698484651e-07, 
2.70599529637464e-17, 5.81535679492809e-17, 2.19440708445348e-15, 
4.34155540862746e-08, -4.06493587341127e-07, -1.61763139817e-06, 
-3.23994151550826e-08, -6.68988064422799e-08, 1.10402768540446e-15, 
-9.91487886840506e-07, -3.94580269988612e-06, -5.8649138992111e-07, 
-5.23882134070119e-06, 1.00998060784975e-16, -1.2029789281118e-07, 
-3.05994985702607e-07, 9.3958523768985e-08, -3.07494963928282e-07, 
-8.46097103856411e-06, -2.86925082960488e-05, 3.08688633134856e-15, 
3.62058033172122e-06, 8.25888178764606e-06, -1.35448644277712e-06, 
-6.34490870510011e-07, -3.72424640639279e-06), Intakerate_kjs = c(-1.48111216166192e-06, 
3.19585611270047e-13, -2.10401730711284e-05, -7.34279384597799e-07, 
1.13218843200315e-16, 2.43314528299791e-16, 9.18139924135334e-15, 
1.81650678296973e-07, -1.70076916943527e-06, -6.76816976994329e-06, 
-1.35559153008866e-07, -2.79904606154499e-07, 4.61925183573226e-15, 
-4.14838531854068e-06, -1.65092384963235e-05, -2.45387997542992e-06, 
-2.19192284894938e-05, 4.22575886324335e-16, -5.03326383521979e-07, 
-1.28028302017971e-06, 3.93122463449433e-07, -1.28655892907593e-06, 
-3.54007028253523e-05, -0.000120049454710668, 1.29155324103624e-14, 
1.51485081079216e-05, 3.45551613995111e-05, -5.66717127657947e-06, 
-2.65470980221389e-06, -1.55822469643474e-05), Intakerate_kjday = c(-0.12796809076759, 
2.76121968137321e-08, -1.81787095334549, -0.0634417388292498, 
9.78210805250721e-12, 2.1022375245102e-11, 7.93272894452929e-10, 
0.0156946186048585, -0.146946456239208, -0.584769868123101, -0.011712310819966, 
-0.0241837579717487, 3.99103358607267e-10, -0.358420491521915, 
-1.42639820608235, -0.212015229877145, -1.89382134149226, 3.65105565784226e-11, 
-0.043487399536299, -0.110616452943527, 0.033965780842031, -0.111158691472161, 
-3.05862072411043, -10.3722728870017, 1.11590200025531e-09, 1.30883110052443, 
2.98556594491776, -0.489643598296466, -0.22936692691128, -1.34630613771962
)), row.names = c(NA, -30L), class = "data.frame")

lat=unique(profit$Latitude)

## for each location I am calculating the weight for    Fall only
nfall=0
latfall<-c(double())
for(i in lat){
  name = paste0("profit.fall.",round(i,5))
  x = subset(profit,Latitude==i & Season=="Fall")
  if(nrow(x)>=1){
    for(j in 1:nrow(x)){
      x$alpha[j]<- 1 # used to be this  x$Density[j]/sum(x$Density)
    }
    nfall= nfall+1
    assign(name, data.frame(x))
    latfall<-c(latfall,round(i,5))
    print(name)
  }
}

View(profit.fall.38.4.959)

profit.fall.all <- data.frame(matrix(ncol=7,nrow=nfall))

names(profit.fall.all)[1]<-'Id'
names(profit.fall.all)[2]<-'Refuge'
names(profit.fall.all)[3]<-'Season'
names(profit.fall.all)[4]<-'HType'
names(profit.fall.all)[5]<-'Lat'
names(profit.fall.all)[6]<-'Long'
names(profit.fall.all)[7]<-'IntakeEnergy'

View(profit.fall.all)

k=0

lat=latfall

for(i in lat){
  df=as.name(paste0('profit.fall.',i))
  d=get(as.character(df))
  x=0
  for(j in 1:nrow(d)){
    x= x + d$Intakerate_kjday[j]*d$alpha[j]
  }
  k=k+1
  new_row <- c(k,as.character(d$Refuge[1]),as.character(d$Season[1]),as.character(d$Habitat.Type[1]),as.numeric(d$Latitude[1]),as.numeric(d$Longtitude[1]),as.numeric(x))
  #names(new_row)<-c("id","Refuge","Season","HType","Lat","Long","Intakerate_kjday")
  #profit.spring.all <- rbind(profit.spring.all, new_row)
  profit.fall.all[k,] <- new_row
}

View(profit.fall.all)

1 Ответ

1 голос
/ 23 апреля 2020

Код, по-видимому, вычисляется (очень неэффективно и неточно) 

sum(d$Intakerate_kjday * d$alpha)

Однако ваша ошибка предполагает, что в одном из фреймов данных отсутствует столбец.

Посмотрите на new_row здесь: 

for(i in lat){
  df=as.name(paste0('profit.fall.',i))
  d=get(as.character(df))
  x=0
  for(j in 1:nrow(d)){
    x= x + d$Intakerate_kjday[j]*d$alpha[j]
  }
  k=k+1
  new_row <- c(k,as.character(d$Refuge[1]),as.character(d$Season[1]),as.character(d$Habitat.Type[1]),as.numeric(d$Latitude[1]),as.numeric(d$Longtitude[1]),as.numeric(x))
  #names(new_row)<-c("id","Refuge","Season","HType","Lat","Long","Intakerate_kjday")
  #profit.spring.all <- rbind(profit.spring.all, new_row)

  if (length(new_row) != ncol(profit.fall.all)) {
    # Catch the bad df
    browser()
  }

  profit.fall.all[k,] <- new_row
}
Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.

Нет похожих вопросов

...