お台場海苔の記録 (2019-2020)
山森光陽（国立教育政策研究所総括研究官（教育心理学））
2020/1/26 更新
0. これまでの記録
1月26日以前の作業内容や海苔の様子については以下を参照
(ファイルサイズの都合で1つのhtmlとしてアップロードできないため)
- 2020年1月3日まで
- 2019年1月16日まで
1. 海苔の成長データ
使うライブラリ
library(gsheet)データの在処は表示しない
- urlにGoogle spreadsheetの編集可リンクを代入する
データの読み込み
data_ <- gsheet2tbl(url)1.1 伸びの記録
1.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)])1.1.2 中間刈り取り後
data.hon <- data_[10:14,]
colnames(data.hon) <- c("Date", "Day", "Temp", "Off_Fuji", "Off_Ctr", "Off_Rbw",
"Bch_Fuji", "Bch_Ctr", "Bch_Rbw")
DT::datatable(data.hon[c(1:9)])1.2 地点別のグラフ
par(family = "HiraKakuProN-W3") #Macintoshの場合
par(oma = c(0, 0, 4, 0))
par(mfrow=c(2,3))
plot(data.md$Day, data.md$Off_Fuji, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "沖・フジTV側")
par(new = T)
plot(data.hon$Day, data.hon$Off_Fuji, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
plot(data.md$Day, data.md$Off_Ctr, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "沖・中央")
par(new = T)
plot(data.hon$Day, data.hon$Off_Ctr, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
plot(data.md$Day, data.md$Off_Rbw, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "沖・レインボー側")
par(new = T)
plot(data.hon$Day, data.hon$Off_Rbw, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
plot(data.md$Day, data.md$Bch_Fuji, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "浜・フジTV側")
par(new = T)
plot(data.hon$Day, data.hon$Bch_Fuji, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
plot(data.md$Day, data.md$Bch_Ctr, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "浜・中央")
par(new = T)
plot(data.hon$Day, data.hon$Bch_Ctr, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
plot(data.md$Day, data.md$Bch_Rbw, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 4, col = 1, xlab = "Day", ylab = "Length (mm) ",
main = "浜・レインボー側")
par(new = T)
plot(data.hon$Day, data.hon$Bch_Rbw, xlim = c(0,50), ylim=c(0,400),
type = "b", pch = 3, col = 1, xlab = "", ylab = "", main = "", axes = F)
2. 海苔の成長予測
2.1 中間刈り取りまでのモデル
\(length = a \times b^{day}\)
2.1.1 中間刈り取りまでの指数回帰
# 指数回帰の係数を推定する
exp.b.f <- nls(Bch_Fuji ~ a*b^Day, data = data.md[1:8,], start = list(a = 1, b = 1))
exp.b.c <- nls(Bch_Ctr ~ a*b^Day, data = data.md[1:8,], start = list(a = 1, b = 1))
exp.b.r <- nls(Bch_Rbw ~ a*b^Day, data = data.md[1:8,], start = list(a = 1, b = 1))
exp.o.f <- nls(Off_Fuji ~ a*b^Day, data = data.md[1:8,], start = list(a = 1, b = 1))
exp.o.c <- nls(Off_Ctr ~ a*b^Day, data = data.md[1:8,], start = list(a = 1, b = 1))
exp.o.r <- nls(Off_Rbw ~ a*b^Day, data = data.md[1:8,], 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.1.2 指数回帰の分析結果
沖・フジTV側
summary(exp.b.f)##
## Formula: Bch_Fuji ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 8.11244 5.73165 1.415 0.207
## b 1.10670 0.02503 44.218 8.96e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 34.61 on 6 degrees of freedom
##
## Number of iterations to convergence: 9
## Achieved convergence tolerance: 2.142e-06沖・中央
summary(exp.b.c)##
## Formula: Bch_Ctr ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 5.96547 3.39517 1.757 0.129
## b 1.11814 0.02015 55.500 2.3e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 25.56 on 6 degrees of freedom
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 5.91e-06沖・レインボー側
summary(exp.b.r)##
## Formula: Bch_Rbw ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 7.34149 5.62367 1.305 0.24
## b 1.10597 0.02714 40.752 1.46e-08 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 33.49 on 6 degrees of freedom
##
## Number of iterations to convergence: 12
## Achieved convergence tolerance: 9.037e-06浜・フジTV側
summary(exp.o.f)##
## Formula: Off_Fuji ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 4.59045 2.39623 1.916 0.104
## b 1.11215 0.01848 60.179 1.41e-09 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 16.06 on 6 degrees of freedom
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 1.068e-06浜・中央
summary(exp.o.c)##
## Formula: Off_Ctr ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 3.21905 1.31137 2.455 0.0495 *
## b 1.12450 0.01443 77.951 3e-10 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 11.2 on 6 degrees of freedom
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 5.497e-07浜・レインボー側
summary(exp.o.r)##
## Formula: Off_Rbw ~ a * b^Day
##
## Parameters:
## Estimate Std. Error t value Pr(>|t|)
## a 2.62538 0.88888 2.954 0.0255 *
## b 1.13659 0.01201 94.622 9.39e-11 ***
## ---
## Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
##
## Residual standard error: 9.7 on 6 degrees of freedom
##
## Number of iterations to convergence: 8
## Achieved convergence tolerance: 6.676e-06例えば，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.2 中間刈り取り後のモデル
\(length = a + day \times b + I({day}^2)\)
2.2.1 多項式近似
# 回帰係数を推定する
af.mlm.b.f <- lm(Bch_Fuji ~ 1 + Day + I(Day^2), data = data.hon)
af.mlm.b.c <- lm(Bch_Ctr ~ 1 + Day + I(Day^2), data = data.hon)
af.mlm.b.r <- lm(Bch_Rbw ~ 1 + Day + I(Day^2), data = data.hon)
af.mlm.o.f <- lm(Off_Fuji ~ 1 + Day + I(Day^2), data = data.hon)
af.mlm.o.c <- lm(Off_Ctr ~ 1 + Day + I(Day^2), data = data.hon)
af.mlm.o.r <- lm(Off_Rbw ~ 1 + Day + I(Day^2), data = data.hon)
# 結果の表示とデータ化
af.sum.mlm.b.f <- summary(af.mlm.b.f); af.sum.mlm.b.f##
## Call:
## lm(formula = Bch_Fuji ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.00e+01 NA NA NA
## Day 2.50e+00 NA NA NA
## I(Day^2) -1.36e-16 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NAaf.sum.mlm.b.c <- summary(af.mlm.b.c); af.sum.mlm.b.c##
## Call:
## lm(formula = Bch_Ctr ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 380.0000 NA NA NA
## Day -21.2500 NA NA NA
## I(Day^2) 0.3125 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NAaf.sum.mlm.b.r <- summary(af.mlm.b.r); af.sum.mlm.b.r##
## Call:
## lm(formula = Bch_Rbw ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 830.000 NA NA NA
## Day -45.000 NA NA NA
## I(Day^2) 0.625 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NAaf.sum.mlm.o.f <- summary(af.mlm.o.f); af.sum.mlm.o.f##
## Call:
## lm(formula = Off_Fuji ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -970.000 NA NA NA
## Day 50.000 NA NA NA
## I(Day^2) -0.625 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NAaf.sum.mlm.o.c <- summary(af.mlm.o.c); af.sum.mlm.o.c##
## Call:
## lm(formula = Off_Ctr ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) -7.00e+01 NA NA NA
## Day 2.50e+00 NA NA NA
## I(Day^2) -1.36e-16 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NAaf.sum.mlm.o.r <- summary(af.mlm.o.r); af.sum.mlm.o.r##
## Call:
## lm(formula = Off_Rbw ~ 1 + Day + I(Day^2), data = data.hon)
##
## Residuals:
## ALL 3 residuals are 0: no residual degrees of freedom!
##
## Coefficients:
## Estimate Std. Error t value Pr(>|t|)
## (Intercept) 1280.0000 NA NA NA
## Day -68.7500 NA NA NA
## I(Day^2) 0.9375 NA NA NA
##
## Residual standard error: NaN on 0 degrees of freedom
## (2 observations deleted due to missingness)
## Multiple R-squared: 1, Adjusted R-squared: NaN
## F-statistic: NaN on 2 and 0 DF, p-value: NA# 推定値の格納
af.b0.b.f <- af.sum.mlm.b.f$coefficients[1,1]; af.b1.b.f <- af.sum.mlm.b.f$coefficients[2,1]; af.b2.b.f <- af.sum.mlm.b.f$coefficients[3,1]
af.b0.b.c <- af.sum.mlm.b.c$coefficients[1,1]; af.b1.b.c <- af.sum.mlm.b.c$coefficients[2,1]; af.b2.b.c <- af.sum.mlm.b.c$coefficients[3,1]
af.b0.b.r <- af.sum.mlm.b.r$coefficients[1,1]; af.b1.b.r <- af.sum.mlm.b.r$coefficients[2,1]; af.b2.b.r <- af.sum.mlm.b.r$coefficients[3,1]
af.b0.o.f <- af.sum.mlm.o.f$coefficients[1,1]; af.b1.o.f <- af.sum.mlm.o.f$coefficients[2,1]; af.b2.o.f <- af.sum.mlm.o.f$coefficients[3,1]
af.b0.o.c <- af.sum.mlm.o.c$coefficients[1,1]; af.b1.o.c <- af.sum.mlm.o.c$coefficients[2,1]; af.b2.o.c <- af.sum.mlm.o.c$coefficients[3,1]
af.b0.o.r <- af.sum.mlm.o.r$coefficients[1,1]; af.b1.o.r <- af.sum.mlm.o.r$coefficients[2,1]; af.b2.o.r <- af.sum.mlm.o.r$coefficients[3,1]2.3 成長予測のグラフ
# 作図する
## 中間刈り取り前 プロット
plot(data.md$Day, data.md$Bch_Fuji, xlim = c(0,50), 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,50), 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,50), 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,50), 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,50), 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,50), ylim = c(0,500),
pch = 2, col = 4, axes = F, xlab = "", ylab = "")
# 中間刈り取り後 プロット
par(new = T)
plot(data.hon$Day, data.hon$Bch_Fuji, xlim = c(0,50), ylim = c(0,500),
pch = 0, col = 2, axes = F, xlab = "", ylab = "")
par(new=T)
plot(data.hon$Day, data.hon$Bch_Ctr, xlim = c(0,50), ylim = c(0,500),
pch = 1, col = 2, axes = F, xlab = "", ylab = "")
par(new=T)
plot(data.hon$Day, data.hon$Bch_Rbw, xlim = c(0,50), ylim = c(0,500),
pch = 2, col = 2, axes = F, xlab = "", ylab = "")
par(new = T)
plot(data.hon$Day, data.hon$Off_Fuji, xlim = c(0,50), ylim = c(0,500),
pch = 0, col = 4, axes = F, xlab="", ylab="")
par(new = T)
plot(data.hon$Day, data.hon$Off_Ctr, xlim = c(0,50), ylim = c(0,500),
pch = 1, col = 4, axes = F, xlab = "", ylab = "")
par(new = T)
plot(data.hon$Day, data.hon$Off_Rbw, xlim = c(0,50), 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)
# 中間刈り取り後 回帰直線
after.day <- c(36:50)
lines(after.day, af.b0.b.f + after.day * af.b1.b.f + af.b2.b.f * (after.day^2), col = 2, lty=1, lwd=2)
lines(after.day, af.b0.b.c + after.day * af.b1.b.c + af.b2.b.c * (after.day^2), col = 2, lty=3, lwd=2)
lines(after.day, af.b0.b.r + after.day * af.b1.b.r + af.b2.b.r * (after.day^2), col = 2, lty=3, lwd=2)
lines(after.day, af.b0.o.f + after.day * af.b1.o.f + af.b2.o.f * (after.day^2), col = 4, lty=1, lwd=2)
lines(after.day, af.b0.o.c + after.day * af.b1.o.c + af.b2.o.c * (after.day^2), col = 4, lty=3, lwd=2)
lines(after.day, af.b0.o.r + after.day * af.b1.o.r + af.b2.o.r * (after.day^2), col = 4, lty=3, 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
)
3 作業の内容
3.1 2019年12月14日
3.2 2019年12月20日
3.3 2019年12月25日
3.5 2020年1月3日
3.6 2020年1月8日
3.7 2020年1月11日
3.9 2020年1月22日
作業内容など
| Jan.22 | |
|---|---|
| 作業時刻 | 21:00-22:00 |
| 水温 | 11 |
| 実測潮位 | 105cm |
| 参加者 | 委員長，ほか5名(台場担当係長不在) |
| 作業内容 | 海水に臭いがあり，濁っている。ほぼ伸びが見られなかったため，網を潮位105cmを基準として，フジテレビ側を-50cm，レインボーブリッジ側を+50cmにして様子を見ることとした。 |
海苔の様子
 |
 |
3.10 2020年1月26日
作業内容など
| Jan.26 | |
|---|---|
| 作業時刻 | 10:00-11:00 |
| 水温 | 10 |
| 実測潮位 | 105cm |
| 参加者 | 委員長，副委員長，小学生1名，ほか6名(台場担当係長不在) |
| 作業内容 | 伸びがあまり見られなかったため，学校側の網を10cm上げた。沖の網はフジテレビ側(105cm-50cm)では育ちが悪く，レインボー側(105cm+50cm)では伸びが良い。検討後，教育用としては良いのではないかとの結論になり，そのままにすることとした。次回1月29日の観察結果次第では，網を上げ干出を多くし白化した部分をなくすようにし，30日に下げることも検討することとした。 |
作業の様子
 |
. |  |
. |  |
海苔の様子
 |
. |  |
. |  |