1 Streszczenie (Executive summary)

Na przestrzeni ostatnich 60 lat da się zauważyć stopniowe karłowacenie śledzi oceanicznych wyławianych w Europie. Analiza danych była procesem wieloetapowym. Z powodu brakujących wartości konieczne okazało się wypełnienie ich przed przejściem do kolejnego etapu. Do ich zastąpienia została wykorzystana średnia wartość ze zbioru pogrupowanych parametrów. W kolejnym etapie przeprowadzona została analiza podstawowych statystyk, a także rozkładu wartości atrybutów i korelacji między nimi. W ostatnim etapie stworzony został regresor przewidujący rozmiar śledzia.

Podczas przeprowadzania analizy zauważono silny wpływ temperatury przy powierzchni wody sst na długość łowionych śledzi length (im wyższa temperatura, tym śledź był mniejszy).

Kolejnymi prametrami, które mogą wpływać na długość są:

  • natężenie połowów fbar,
  • łączna liczba ryb złowionych w ramach połowu totaln,
  • roczny narybek recr.

Przyczyną tego mogą być działania ludzi, czyli zwiększenie liczby wyławianych śledzi, co spowodowało opóźnienie ich rozwoju.

2 Wstęp

Celem projektu jest określenie jakie mogą być główne przyczyny stopniowego karłowacenia śledzi oceanicznych wyławianych w Europie. Dane do analizy zawierają pomiary śledzi i warunków w jakich żyją z ostatnich 60 lat, były one pobierane z połowów komercyjnych jednostek. W ramach połowu jednej jednostki losowo wybierano od 50 do 100 sztuk trzyletnich śledzi. Dane są uporządkowane chronologicznie.

Kolejne kolumny w zbiorze danych to:

  • length: długość złowionego śledzia [cm];
  • cfin1: dostępność planktonu [zagęszczenie Calanus finmarchicus gat. 1];
  • cfin2: dostępność planktonu [zagęszczenie Calanus finmarchicus gat. 2];
  • chel1: dostępność planktonu [zagęszczenie Calanus helgolandicus gat. 1];
  • chel2: dostępność planktonu [zagęszczenie Calanus helgolandicus gat. 2];
  • lcop1: dostępność planktonu [zagęszczenie widłonogów gat. 1];
  • lcop2: dostępność planktonu [zagęszczenie widłonogów gat. 2];
  • fbar: natężenie połowów w regionie [ułamek pozostawionego narybku];
  • recr: roczny narybek [liczba śledzi];
  • cumf: łączne roczne natężenie połowów w regionie [ułamek pozostawionego narybku];
  • totaln: łączna liczba ryb złowionych w ramach połowu [liczba śledzi];
  • sst: temperatura przy powierzchni wody [°C];
  • sal: poziom zasolenia wody [Knudsen ppt];
  • xmonth: miesiąc połowu [numer miesiąca];
  • nao: oscylacja północnoatlantycka [mb].

3 Wstępne przetwarzanie

W tej części następuje załadowanie wykorzystanych w projekcie bibliotek, zapewnienie powtarzalności wyników, wczytywanie danych i wstępne przygotowanie danych.

3.1 Wykorzystanie biblioteki

library(dplyr) 
library(tidyr) 
library(skimr) 
library(ggplot2)
library(plotly) 
library(ggcorrplot) 
library(caret)

3.2 Powtarzalność wyników

set.seed(23)

3.3 Wczytywanie danych

Wartości puste w zbiorze danych oznaczone są znakiem ‘?’, przy wczytywaniu, korzystając z parametru na.string, zostaną zastąpione wartością NA.

df <- read.csv('sledzie.csv', na.string='?')
head(df)
dim(df)
## [1] 52582    16

Zbiór danych zawiera 52582 obserwacji opisanych 16 cechami.

3.4 Przetwarzanie brakujących danych

Z powyższego wykresu widać, że brakujące wartości występują w kolumnach, które opisują:

  • dostępność planktonu (cfin1, cfin2, chel1, chel2, lcop1, lcop2),
  • temperaturę przy powierzchni wody (sst).

W przypadku odrzucenia tych pomiarów utracone zostanie zbyt dużo istotnych informacji. Można zauważyć, że dane te są silnie ze sobą powiązane (dla danego połowu przyjmują stałe wartości).

df[200:215,]

Jak widać w powyższym zestawieniu, brakujące dane można zastąpić danymi z podobnych wierszy (obserwacji).

impute.mean <- function(x) replace(x, is.na(x), mean(x, na.rm = TRUE))

df <- df %>%
          group_by(recr, xmonth) %>%
          mutate(
             cfin1 = impute.mean(cfin1),
             cfin2 = impute.mean(cfin2),
             chel1 = impute.mean(chel1),
             chel2 = impute.mean(chel2),
             lcop1 = impute.mean(lcop1),
             lcop2 = impute.mean(lcop2),
             sst = impute.mean(sst) 
          ) %>% ungroup()

4 Opis danych

W tej części następuje opisane głównych statystyk zbioru, dokonanie analizy wartości atrybutów oraz przedstawienie wykresu zagęszczenia planktonu w czasie.

Wartość
Liczba kolumn (zmiennych) 16
Liczba wierszy (obserwacji) 52582 
prettyTable(skim(df), show_entries=16)

4.1 Analiza wartości atrybutów

length

Na poniższym wykresie można z łatwością zauważyć, że rozkład długości złowionego śledzia length przypomina rozkład normalny.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   19.00   24.00   25.50   25.30   26.50   32.50    1.65

cfin1

Z wykresu wynika, że dostępność planktonu cfin1 była zdecydowanie niska.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##  0.0000  0.0000  0.1111  0.4457  0.3333 37.6667  0.9500

cfin2

Podobnie jak w przypadku cfin1, rozkład cfin2 jest zdominowany przez niskie wartości, jednak w znacznie mniejszym stopniu.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##  0.0000  0.2778  0.7012  2.0258  1.7936 19.3958  3.7100

chel1

Podobnie jak w poprzednich cechach można zauważyć niską dostępność planktonu chel1.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   0.000   2.469   5.750  10.003  11.500  75.000  14.270

chel2

Można zauważyć, że rozkład dostępności planktonu chel2 jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   5.238  13.427  21.435  21.219  27.193  57.706   9.980

lcop1

Można zauważyć, że rozkład dostępności planktonu lcop1 jest zróżnicowany.

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max.       SD 
##   0.3074   2.5479   7.0717  12.8048  21.2315 115.5833  15.0400

lcop2

Widać, że rozkład dostępności planktonu lcop2 jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   7.849  17.808  24.859  28.422  37.232  68.736  13.850

fbar

Widać, że rozkład natężenia połowów w regionie fbar jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##  0.0680  0.2270  0.3320  0.3304  0.4560  0.8490  0.1600

recr

Widać, że rozkład rocznego narybeku recr jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##  140515  360061  421391  520367  724151 1565890  271063

cumf

Widać, że rozkład łącznego roczne natężenia połowów w regionie cumf jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
## 0.06833 0.14809 0.23191 0.22981 0.29803 0.39801 0.09000

totaln

Widać, że rozkład łącznej liczby ryb złowionych w ramach połowu totaln jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##  144137  306068  539558  514973  730351 1015595  221382

sst

Widać, że rozkład temperatury przy powierzchni wody sst jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   12.77   13.60   13.86   13.87   14.16   14.73    0.44

sal

Widać, że rozkład poziomu zasolenia wody sal jest zróżnicowany.

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   35.40   35.51   35.51   35.51   35.52   35.61    0.04

xmonth

Z wykresu można odczytać, że większość obserwacji pochodzi z okresu letnio-jesiennego (lipiec, sierpień, wrzesień, październik).

##    Min. 1st Qu.  Median    Mean 3rd Qu.    Max.      SD 
##   1.000   5.000   8.000   7.258   9.000  12.000   2.760

nao

Widać, że rozkład oscylacji północnoatlantyckiej nao jest zróżnicowany.

##     Min.  1st Qu.   Median     Mean  3rd Qu.     Max.       SD 
## -4.89000 -1.89000  0.20000 -0.09236  1.63000  5.08000  2.25000

4.2 Analiza dostępności planktonów w czasie

Wykres przedstawia dostępność planktonów w czasie.

5 Analiza korelacji między zmiennymi

W tej części następuje analiza korelacji między zmiennymi.

corrMatrix <- round( cor(df %>% select(-X)), 2)
ggcorrplot(round(corrMatrix,1), type = "lower", lab = TRUE)

Bardzo silne korelacje dodatnie:

  • 1.0 między chel1, a lcop1,
  • 0.9 między chel2, a lcop2.

Silne korelacje dodatnie:

  • 0.8 między fbar, a cumf (oba parametry opisują natężenie połowów w regionie)
  • 0.7 między cfin2, a lcop2

Silna korelacja ujemna:

  • -0.7 między totaln, a cumf (parametry opisują liczbę złowionych śledzi i ułamek pozostawionego narybku)

Do oszacowania długości śledzia length istotny wydaje się być parametr temperatury przy powierzchni wody sst, oraz w mniejszym stopniu parametry: oscylacja północnoatlantycka nao, zagęszczenie Calanus helgolandicus gat. 1 chel1, zagęszczenie widłonogów gat. 1 lcop1 oraz ułamek pozostawionego narybku fbar. Zmienne te warto brać pod uwagę podczas analizy karłowacenia się śledzi.

5.1 Korelacja:

length~sst (-0.45)

length~nao (-0.26)

length~chel1 (0.22)

length~lcop1 (0.24)

length~fbar (0.25)

chel1~lcop1 (0.95)

chel2~lcop2 (0.89)

fbar~cumf (0.82)

cfin2~lcop2 (0.65)

totaln~cumf (-0.71)

6 Zmiana rozmiaru śledzi w czasie

Wykres przedstawia zmianę rozmiaru złowionych śledzi na przestrzeni lat.

7 Przewidywanie rozmiaru śledzia

W tej części zostanie zaproponowany przykładowy model do predykcji długości śledzia. Zbiór został podzielony na dane uczące, walidujące i testowe w proporcjach 14:3:3.

7.1 Wyniki predykcji

Do modelu został wykorzystany algorytm Random Forest, a najlepsze wyniki uzyskano dla modelu o parametrze mtry = 4.

inTrainging <- createDataPartition(df$length,
                                   p = .7,
                                   list = FALSE)

training <- df[ inTrainging, ] %>% select(-c(X, xmonth))
testing  <- df[-inTrainging, ] 

ctrl <- trainControl(
    method = "repeatedcv",
    number = 2,
    repeats = 5
    )

fit <- train(length ~ .,
             data = training,
             method = "rf",
             trControl = ctrl,
             tuneGrid = expand.grid(mtry = 4:4), 
             ntree = 32
             )
fit
## Random Forest 
## 
## 36810 samples
##    13 predictor
## 
## No pre-processing
## Resampling: Cross-Validated (2 fold, repeated 5 times) 
## Summary of sample sizes: 18405, 18405, 18405, 18405, 18405, 18405, ... 
## Resampling results:
## 
##   RMSE      Rsquared   MAE      
##   1.187991  0.4840271  0.9401704
## 
## Tuning parameter 'mtry' was held constant at a value of 4
knitr::kable(fit$results[,1:3], caption="Tabela wyników ze zbioru uczącego")
Tabela wyników ze zbioru uczącego
mtry RMSE Rsquared
4 1.187991 0.4840271 
predictions <- predict(fit, newdata = testing %>% select(-c(X, xmonth)))
result_predictions <- postResample(pred = predictions, obs=testing$length)
knitr::kable(result_predictions[1:2], caption = "Tabela wyników ze zbioru testującego")
Tabela wyników ze zbioru testującego
x
RMSE 1.1873360
Rsquared 0.4826705

7.2 Wykres przedstawiający porównanie wartości przewidzianych z rzeczywistymi

7.3 Analizę ważności atrybutów najlepszego znalezionego modelu regresji.

W tej części zostanie przedstawiona analiza ważności atrybutów modelu regresji.

ggplotly( 
  ggplot(varImp(fit)) + geom_col(fill = 'steelblue') + theme_minimal() + 
  labs(y="Ważność", x="Cecha") )

Z wykresu wynika, że parametry: temperatura przy powierzchni wody sst, natężenie połowów fbar, łączna liczba ryb złowionych w ramach połowu totaln oraz roczny narybek recr zdecydowanie wyróżniają się od pozostałych parametrów. Zmiany temperatury mogą powodować zmiany w rozwoju planktonów, które mogą pogarszać warunki środowiska życia śledzi i ograniczać ich wzrost. Logiczne wydaje się, że jeśli zwiększymy liczbę wyławianych śledzi to z biegiem czasu doprowadzi do zmniejszenia ich populacji, co wiąże się ze zmniejszeniem ich rozwoju do czasu kolejnego połowu.