1. 海苔の成長データ

使うライブラリ

library(gsheet)

データの在処は表示しない

データの読み込み

data_ <- gsheet2tbl(url)

1.1 中間刈り取りまでのデータ

data.md <- data_[1:9,]
colnames(data.md) <- c("Date", "Day", "Temp", "Off_Fuji", "Off_Ctr", "Off_Rbw",
                       "Bch_Fuji", "Bch_Ctr", "Bch_Rbw")
DT::datatable(data.md[c(1:9)])

表の凡例 Temp: Water temperature, Off: Offshore, Bch: Beach, Fuji: Fuji TV, Ctr: Center, Rbw: Rainbow bridge

1.2 地点別のグラフ

par(family = "HiraKakuProN-W3") #Macintoshの場合
par(oma = c(0, 0, 4, 0)) 
par(mfrow=c(2,3)) 

plot(data_$Day, data_$Offshore_Fuji, xlim = c(0,50), ylim=c(0,500), 
     type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ", 
     main = "沖・フジTV側")

plot(data_$Day, data_$Offshore_Center, xlim = c(0,50), ylim = c(0,500), 
     type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
     main = "沖・中央")

plot(data_$Day, data_$Offshore_Rainbow, xlim = c(0,50), ylim = c(0,500), 
     type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
     main = "沖・レインボー側")

plot(data_$Day, data_$Beach_Fuji, xlim=c(0,50), ylim=c(0,500), type = "b", 
     pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ", 
     main = "浜・フジTV側")

plot(data_$Day, data_$Beach_Center, xlim=c(0,50), ylim=c(0,500), 
     type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ", 
     main = "浜・中央")

plot(data_$Day, data_$Beach_Rainbow, xlim = c(0,50), ylim = c(0,500), 
     type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ", 
     main = "浜・レインボー側")

2. 海苔の成長予測

2.1 中間刈り取りまでのモデル

\(length = a \times b^{day}\)

2.2 中間刈り取りまでの指数回帰

# 指数回帰の係数を推定する
exp.b.f <- nls(Bch_Fuji ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 
exp.b.c <- nls(Bch_Ctr  ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 
exp.b.r <- nls(Bch_Rbw  ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 

exp.o.f <- nls(Off_Fuji ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 
exp.o.c <- nls(Off_Ctr  ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 
exp.o.r <- nls(Off_Rbw  ~ a*b^Day, data = data.md, start = list(a = 1, b = 1)) 

# 推定値を格納する
sum.b.f <- summary(exp.b.f)
sum.b.c <- summary(exp.b.c)
sum.b.r <- summary(exp.b.r)

sum.o.f <- summary(exp.o.f)
sum.o.c <- summary(exp.o.c)
sum.o.r <- summary(exp.o.r)

b.f.a <- sum.b.f$coefficients[1,1]; b.f.b <- sum.b.f$coefficients[2,1]
b.c.a <- sum.b.c$coefficients[1,1]; b.c.b <- sum.b.c$coefficients[2,1]
b.r.a <- sum.b.r$coefficients[1,1]; b.r.b <- sum.b.r$coefficients[2,1]

o.f.a <- sum.o.f$coefficients[1,1]; o.f.b <- sum.o.f$coefficients[2,1]
o.c.a <- sum.o.c$coefficients[1,1]; o.c.b <- sum.o.c$coefficients[2,1]
o.r.a <- sum.o.r$coefficients[1,1]; o.r.b <- sum.o.r$coefficients[2,1]

2.3 指数回帰の分析結果

沖・フジTV側

summary(exp.b.f)
## 
## Formula: Bch_Fuji ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a   0.5109     0.5133   0.995   0.5015  
## b   1.1632     0.1053  11.050   0.0575 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8063 on 1 degrees of freedom
## 
## Number of iterations to convergence: 11 
## Achieved convergence tolerance: 3.574e-06
##   (6 observations deleted due to missingness)

沖・中央

summary(exp.b.c)
## 
## Formula: Bch_Ctr ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a  0.37882    0.25379   1.493   0.3758  
## b  1.24096    0.07205  17.223   0.0369 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5513 on 1 degrees of freedom
## 
## Number of iterations to convergence: 7 
## Achieved convergence tolerance: 4.653e-06
##   (6 observations deleted due to missingness)

沖・レインボー側

summary(exp.b.r)
## 
## Formula: Bch_Rbw ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a  0.37882    0.25379   1.493   0.3758  
## b  1.24096    0.07205  17.223   0.0369 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5513 on 1 degrees of freedom
## 
## Number of iterations to convergence: 7 
## Achieved convergence tolerance: 4.653e-06
##   (6 observations deleted due to missingness)

浜・フジTV側

summary(exp.o.f)
## 
## Formula: Off_Fuji ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a   0.5109     0.5133   0.995   0.5015  
## b   1.1632     0.1053  11.050   0.0575 .
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.8063 on 1 degrees of freedom
## 
## Number of iterations to convergence: 11 
## Achieved convergence tolerance: 3.574e-06
##   (6 observations deleted due to missingness)

浜・中央

summary(exp.o.c)
## 
## Formula: Off_Ctr ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a  0.37882    0.25379   1.493   0.3758  
## b  1.24096    0.07205  17.223   0.0369 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5513 on 1 degrees of freedom
## 
## Number of iterations to convergence: 7 
## Achieved convergence tolerance: 4.653e-06
##   (6 observations deleted due to missingness)

浜・レインボー側

summary(exp.o.r)
## 
## Formula: Off_Rbw ~ a * b^Day
## 
## Parameters:
##   Estimate Std. Error t value Pr(>|t|)  
## a  0.37882    0.25379   1.493   0.3758  
## b  1.24096    0.07205  17.223   0.0369 *
## ---
## Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
## 
## Residual standard error: 0.5513 on 1 degrees of freedom
## 
## Number of iterations to convergence: 7 
## Achieved convergence tolerance: 4.653e-06
##   (6 observations deleted due to missingness)

これらの結果をまとめると以下の通り

海苔の長さの推定値を求める方法

沖・フジTV側

沖・中央
沖・レインボー側
浜・フジTV側
浜・中央
浜・レインボー側

例えば，20日後の長さの推定値を求めたい場合，\(length = a \times b^{day}\)であるから，aの値と，bの値の20乗をかけることで，長さの推定値が求まる。

仮に\(a = 0.38\), \(b = 1.24\)で20日後の長さを求めるのであれば，\(0.38 \times 1.24^{20} = 28.06\)となる。

2.4 成長予測のグラフ

# 作図する
plot(data.md$Day, data.md$Bch_Fuji, xlim = c(0,36), ylim = c(0,500), 
     pch = 0, col = 2,xlab = "Day", ylab = "Length (mm) ",
     main = "Odaiba-nori Growth Estimate 2019-20")
par(new=T)
plot(data.md$Day, data.md$Bch_Ctr, xlim = c(0,36), ylim = c(0,500),
     pch = 1, col = 2, axes = F, xlab = "", ylab = "")
par(new=T)
plot(data.md$Day, data.md$Bch_Rbw, xlim = c(0,36), ylim = c(0,500), 
     pch = 2, col = 2, axes = F, xlab = "", ylab = "")
par(new = T)

plot(data.md$Day, data.md$Off_Fuji, xlim = c(0,36), ylim = c(0,500), 
     pch = 0, col = 4, axes = F, xlab="", ylab="")
par(new = T)
plot(data.md$Day, data.md$Off_Ctr, xlim = c(0,36), ylim = c(0,500),
     pch = 1, col = 4, axes = F, xlab = "", ylab = "")
par(new = T)
plot(data.md$Day, data.md$Off_Rbw, xlim = c(0,36), ylim = c(0,500), 
     pch = 2, col = 4, axes = F, xlab = "", ylab = "")

data.day <- c(1:36)
lines(data.day, b.f.a*b.f.b^data.day, col = 2, lty = 1, lwd = 2)
lines(data.day, b.c.a*b.c.b^data.day, col = 2, lty = 3, lwd = 2)
lines(data.day, b.r.a*b.r.b^data.day, col = 2, lty = 4, lwd = 2)

lines(data.day, o.f.a*o.f.b^data.day, col = 4, lty = 1, lwd = 2)
lines(data.day, o.c.a*o.c.b^data.day, col = 4, lty = 3, lwd = 2)
lines(data.day, o.r.a*o.r.b^data.day, col = 4, lty = 4, lwd = 2)

legend (0, 500, 
c("Beach_Fuji", "Beach_Center", "Beach_Rainbow", 
  "Offshore_Fuji", "Offshore_Center", "Offshore_Rainbow"),
lwd = 1, 
pch = c(0,1,2,0,1,2), lty = c(1,3,4,1,3,4), col = c(2,2,2,4,4,4), 
bty = "n", bg = "n", cex = 0.8
)

  • このグラフは凡例の順に描画しているため回帰直線が重なり見えないことがある

2.5 観察・作業予定日ごとの成長予測

# dec 29 (16days)
days <- 16
d29.b.f <- b.f.a*b.f.b^days; d29.b.c <- b.c.a*b.c.b^days
d29.b.r <- b.r.a*b.r.b^days; d29.o.f <- o.f.a*o.f.b^days
d29.o.c <- o.c.a*o.c.b^days; d29.o.r <- o.r.a*o.r.b^days

# Jan 3 (21days)
days <- 21
j3.b.f <- b.f.a*b.f.b^days; j3.b.c <- b.c.a*b.c.b^days
j3.b.r <- b.r.a*b.r.b^days; j3.o.f <- o.f.a*o.f.b^days
j3.o.c <- o.c.a*o.c.b^days; j3.o.r <- o.r.a*o.r.b^days

# Jan 8 (26days)
days <- 26
j8.b.f <- b.f.a*b.f.b^days; j8.b.c <- b.c.a*b.c.b^days
j8.b.r <- b.r.a*b.r.b^days; j8.o.f <- o.f.a*o.f.b^days
j8.o.c <- o.c.a*o.c.b^days; j8.o.r <- o.r.a*o.r.b^days

# Jan 11 (29days)
days <- 29
j11.b.f <- b.f.a*b.f.b^days; j11.b.c <- b.c.a*b.c.b^days
j11.b.r <- b.r.a*b.r.b^days; j11.o.f <- o.f.a*o.f.b^days
j11.o.c <- o.c.a*o.c.b^days; j11.o.r <- o.r.a*o.r.b^days

# Jan 16 (34days)
days <- 34
j16.b.f <- b.f.a*b.f.b^days; j16.b.c <- b.c.a*b.c.b^days
j16.b.r <- b.r.a*b.r.b^days; j16.o.f <- o.f.a*o.f.b^days
j16.o.c <- o.c.a*o.c.b^days; j16.o.r <- o.r.a*o.r.b^days

# Jan 18 (36days)
days <- 36
j18.b.f <- b.f.a*b.f.b^days; j18.b.c <- b.c.a*b.c.b^days
j18.b.r <- b.r.a*b.r.b^days; j18.o.f <- o.f.a*o.f.b^days
j18.o.c <- o.c.a*o.c.b^days; j18.o.r <- o.r.a*o.r.b^days

md.es.tab <- data.frame(t(
             matrix(c(d29.b.f, d29.b.c, d29.b.r, d29.o.f, d29.o.c, d29.o.r,
                        j3.b.f,  j3.b.c,  j3.b.r,  j3.o.f,  j3.o.c,  j3.o.r,
                        j8.b.f,  j8.b.c,  j8.b.r,  j8.o.f,  j8.o.c,  j8.o.r,
                       j11.b.f, j11.b.c, j11.b.r, j11.o.f, j11.o.c, j11.o.r,
                       j16.b.f, j16.b.c, j16.b.r, j16.o.f, j16.o.c, j16.o.r,
                       j18.b.f, j18.b.c, j18.b.r, j18.o.f, j18.o.c, j18.o.r),
                    nrow = 6, ncol = 6)
             ))

colnames(md.es.tab) <- c( "Off_Fuji", "Off_Ctr", "Off_Rbw",
                          "Bch_Fuji", "Bch_Ctr", "Bch_Rbw")
rownames(md.es.tab) <- c("Dec.29", "Jan.03", "Jan.08", 
                         "Jan.11", "Jan.16", "Jan.18")

DT::datatable(format(md.es.tab, digits = 2, nsmall = 2))
  • 単位はmm

3 作業の内容

3.1 2019年12月14日

作業内容など

Dec.14
作業時刻 22:00-22-30
水温 NA
実測潮位 47cm
参加者 NA
作業内容 沖に張った海苔網を1週間海面上にならないように下げた。

3.2 2019年12月20日

作業内容など

Dec.20
作業時刻 19:00-12:10
水温 NA
実測潮位 107cm
参加者 委員長，副委員長，台場担当係長ほか6名
作業内容 沖に張った海苔網6枚のうち3枚を学校側のひびに移設。沖，学校側ともに潮位105cmの高さに張り直し。

作業の様子

. .

海苔の様子

3.3 2019年12月25日

作業内容など

Dec.25
作業時刻 19:00-19:40
水温 14
実測潮位 108cm
参加者 委員長，校長，副校長，5年担任，4年担任ほか9名(台場担当係長不在)
作業内容 観察。学校側・フジテレビ側の網がやや高かったため海面より3cm下げた(105cmにそろえた)。

作業の様子

.

海苔の様子

. .