Grau de Ciències del Mar · Universitat de Barcelona
Autor/a
Pràctiques d’Ecologia Marina
Publicat
14 de maig de 2026
Com fer servir aquest document
Executa cada bloc de codi en ordre, de dalt a baix. Llegeix les explicacions de cada secció abans d’executar el codi. Tots els gràfics es guardaran automàticament a la carpeta grafics_visual/.
1 Càrrega i preparació de dades
1.1 Importació del fitxer CSV
El primer pas és sempre carregar les dades i revisar-ne l’estructura. Usem read.csv2() perquè el fitxer utilitza ; com a separador (format europeu) i , com a separador decimal.
Mostra/amaga el codi
# Llegim el CSV i netegem l'encodingdades_raw <-read.csv2("ECONSHAB_Peixos_2026.csv",stringsAsFactors =FALSE,fileEncoding ="latin1") %>%mutate(across(where(is.character), \(x) iconv(x, from ="latin1", to ="UTF-8")))# Comprovem les primeres files i l'estructuraglimpse(dades_raw)
L’exploració inicial ens permet entendre l’estructura del dataset, identificar valors mancants i detectar possibles outliers abans de fer cap anàlisi.
2.1 Resum estadístic
Mostra/amaga el codi
# Resum de les variables numèriques principalsresum <- dades %>%select(nombre, biomassa, talla.mediana, prof, pendent, rugositat) %>%summary()print(resum)
nombre biomassa talla.mediana prof
Min. : 1.000 Min. : 1.26 Min. : 4.00 Min. : 2.0
1st Qu.: 1.000 1st Qu.: 15.51 1st Qu.: 7.00 1st Qu.: 2.0
Median : 3.000 Median : 149.80 Median :15.00 Median : 2.5
Mean : 8.828 Mean : 850.13 Mean :16.93 Mean : 1319.7
3rd Qu.: 6.000 3rd Qu.: 645.87 3rd Qu.:22.00 3rd Qu.: 3.0
Max. :300.000 Max. :34058.13 Max. :70.00 Max. :46115.0
pendent rugositat
Min. :1.000 Min. : 1.5
1st Qu.:1.000 1st Qu.: 2.0
Median :1.500 Median : 2.0
Mean :1.517 Mean : 798.4
3rd Qu.:1.500 3rd Qu.: 2.5
Max. :3.000 Max. :46083.0
Mostra/amaga el codi
# Resum per grup de protecciódades %>%mutate(densitat = nombre /250) %>%group_by(proteccio) %>%summarise(n_obs =n(),dens_mitja =round(mean(densitat, na.rm =TRUE), 4),dens_sd =round(sd(densitat, na.rm =TRUE), 4),bio_mitja =round(mean(biomassa, na.rm =TRUE), 2),bio_sd =round(sd(biomassa, na.rm =TRUE), 2),n_especies =n_distinct(especie),.groups ="drop" ) %>%kable(caption ="Resum per grau de protecció (Visual)",col.names =c("Protecció", "N obs.", "Dens. mitja", "Dens. SD","Bio. mitja (g)", "Bio. SD", "N espècies")) %>%kable_styling(bootstrap_options =c("striped", "hover", "condensed"))
Resum per grau de protecció (Visual)
Protecció
N obs.
Dens. mitja
Dens. SD
Bio. mitja (g)
Bio. SD
N espècies
AMP
1159
0.0336
0.0878
919.37
2450.99
21
NP
346
0.0412
0.0896
618.19
2604.51
17
2.2 Boxplot densitat AMP i NP
Consell
Què veiem? El boxplot per espècie permet comparar la distribució de la densitat dins de cada grup de protecció. Les caixes representen el rang interquartílic (IQR), la línia central la mediana, i els punts els outliers.
Mostra/amaga el codi
# Preparem les dades de densitatdades_dens <- dades %>%mutate(densitat = nombre /250)# Boxplot per AMPp_box_dens_AMP <- dades_dens %>%filter(proteccio =="AMP") %>%ggplot(aes(x =reorder(especie, densitat, FUN = median),y = densitat, fill = especie)) +geom_boxplot(outlier.shape =21, outlier.size =1.5, show.legend =FALSE) +labs(title ="Densitat per espècie – AMP (Àrea Marina Protegida)",subtitle ="Tècnica: Comptatge visual | Unitat: ind/250m²",x =NULL, y ="Densitat (ind/250m²)") +theme_bw(base_size =12) +theme(axis.text.x =element_text(angle =45, hjust =1, face ="italic"),panel.grid.major.x =element_blank())print(p_box_dens_AMP)
Mostra/amaga el codi
desa_grafic("01a_boxplot_densitat_AMP.png", p_box_dens_AMP, ample =12, alt =7)
Mostra/amaga el codi
# Boxplot per NPp_box_dens_NP <- dades_dens %>%filter(proteccio =="NP") %>%ggplot(aes(x =reorder(especie, densitat, FUN = median),y = densitat, fill = especie)) +geom_boxplot(outlier.shape =21, outlier.size =1.5, show.legend =FALSE) +labs(title ="Densitat per espècie – NP (No Protegit)",subtitle ="Tècnica: Comptatge visual | Unitat: ind/250m²",x =NULL, y ="Densitat (ind/250m²)") +theme_bw(base_size =12) +theme(axis.text.x =element_text(angle =45, hjust =1, face ="italic"),panel.grid.major.x =element_blank())print(p_box_dens_NP)
Mostra/amaga el codi
desa_grafic("01b_boxplot_densitat_NP.png", p_box_dens_NP, ample =12, alt =7)
2.3 Boxplot biomassa AMP i NP
Mostra/amaga el codi
dades_bio <- dades %>%mutate(biomassa_250 = biomassa /250)p_box_bio_AMP <- dades_bio %>%filter(proteccio =="AMP") %>%ggplot(aes(x =reorder(especie, biomassa_250, FUN = median),y = biomassa_250, fill = especie)) +geom_boxplot(outlier.shape =21, outlier.size =1.5, show.legend =FALSE) +labs(title ="Biomassa per espècie – AMP",x =NULL, y ="Biomassa (g/250m²)") +theme_bw(base_size =12) +theme(axis.text.x =element_text(angle =45, hjust =1, face ="italic"),panel.grid.major.x =element_blank())p_box_bio_NP <- dades_bio %>%filter(proteccio =="NP") %>%ggplot(aes(x =reorder(especie, biomassa_250, FUN = median),y = biomassa_250, fill = especie)) +geom_boxplot(outlier.shape =21, outlier.size =1.5, show.legend =FALSE) +labs(title ="Biomassa per espècie – NP",x =NULL, y ="Biomassa (g/250m²)") +theme_bw(base_size =12) +theme(axis.text.x =element_text(angle =45, hjust =1, face ="italic"),panel.grid.major.x =element_blank())print(p_box_bio_AMP)
Mostra/amaga el codi
print(p_box_bio_NP)
Mostra/amaga el codi
desa_grafic("01c_boxplot_biomassa_AMP.png", p_box_bio_AMP, ample =12, alt =7)desa_grafic("01d_boxplot_biomassa_NP.png", p_box_bio_NP, ample =12, alt =7)
3 Diversitat
Conceptes clau
S (riquesa): nombre d’espècies observades per transsecte.
H’ (Shannon): mesura la diversitat tenint en compte la riquesa i l’equitativitat. Valors alts → comunitat diversa i equilibrada.
J’ (Pielou): equitativitat normalitzada entre 0 i 1. J’ = H’ / ln(S). Valors propers a 1 → les espècies estan distribuïdes uniformement.
3.1 Preparació de la matriu d’espècies × transsecte
p_riquesa <-ggplot(riquesa_df, aes(x = proteccio, y = S, fill = proteccio)) +geom_boxplot(alpha =0.7, outlier.shape =21) +geom_jitter(width =0.1, size =2.5, alpha =0.7) +scale_fill_manual(values = colors_prot) +labs(title ="Riquesa d'espècies (S) per transsecte",subtitle ="Cada punt = un transsecte",x =NULL, y ="Nombre d'espècies (S)") +theme_bw(base_size =13) +theme(legend.position ="none")print(p_riquesa)
Mostra/amaga el codi
desa_grafic("02a_boxplot_riquesa_S.png", p_riquesa, ample =6, alt =5)
3.3 Índex de Shannon (H’) i equitativitat de Pielou (J’)
Mostra/amaga el codi
H <-diversity(mat_num, index ="shannon")S <-specnumber(mat_num)J <- H /log(S)div_df <-tibble(transsecte =rownames(mat_num),H_shannon = H,S = S,J_pielou = J) %>%left_join(trans_sp %>%distinct(transsecte, proteccio, estacio), by ="transsecte")# Taula resumdiv_df %>%group_by(proteccio) %>%summarise(H_mitja =round(mean(H_shannon, na.rm =TRUE), 3),H_sd =round(sd(H_shannon, na.rm =TRUE), 3),J_mitja =round(mean(J_pielou, na.rm =TRUE), 3),J_sd =round(sd(J_pielou, na.rm =TRUE), 3),.groups ="drop" ) %>%kable(caption ="Índexs de diversitat per grau de protecció",col.names =c("Protecció", "H' mitja", "H' SD", "J' mitja", "J' SD")) %>%kable_styling(bootstrap_options =c("striped", "hover"))
Índexs de diversitat per grau de protecció
Protecció
H' mitja
H' SD
J' mitja
J' SD
AMP
1.431
0.432
0.689
0.19
NP
1.299
0.312
0.713
0.17
Mostra/amaga el codi
p_shannon <-ggplot(div_df, aes(x = proteccio, y = H_shannon, fill = proteccio)) +geom_boxplot(alpha =0.7, outlier.shape =21) +geom_jitter(width =0.1, size =2.5, alpha =0.7) +scale_fill_manual(values = colors_prot) +labs(title ="Índex de Shannon (H')", x =NULL, y ="H'") +theme_bw(base_size =13) +theme(legend.position ="none")p_pielou <-ggplot(div_df, aes(x = proteccio, y = J_pielou, fill = proteccio)) +geom_boxplot(alpha =0.7, outlier.shape =21) +geom_jitter(width =0.1, size =2.5, alpha =0.7) +scale_fill_manual(values = colors_prot) +labs(title ="Equitativitat de Pielou (J')", x =NULL, y ="J'") +theme_bw(base_size =13) +theme(legend.position ="none")p_shannon + p_pielou
Mostra/amaga el codi
desa_grafic("02b_boxplot_shannon.png", p_shannon, ample =5, alt =5)desa_grafic("02c_boxplot_pielou.png", p_pielou, ample =5, alt =5)
3.4 Tests estadístics de comparació
Consell
Per quin test optem? Primer comprovem la normalitat (Shapiro-Wilk). Si p < 0.05 (dades no normals), usem el test de Wilcoxon (2 grups) o Kruskal-Wallis (> 2 grups). Si les dades són normals, podríem usar una t de Student o ANOVA.
-- S --
AMP: n=87, mediana=8.000
NP: n=40, mediana=6.000
Wilcoxon rank sum test with continuity correction
data: amp_vals and np_vals
W = 2718.5, p-value = 2.689e-07
alternative hypothesis: true location shift is not equal to 0
-- H_shannon --
AMP: n=87, mediana=1.539
NP: n=40, mediana=1.344
Wilcoxon rank sum test with continuity correction
data: amp_vals and np_vals
W = 2240, p-value = 0.009526
alternative hypothesis: true location shift is not equal to 0
-- J_pielou --
AMP: n=87, mediana=0.719
NP: n=40, mediana=0.757
Wilcoxon rank sum test with continuity correction
data: amp_vals and np_vals
W = 1630, p-value = 0.5698
alternative hypothesis: true location shift is not equal to 0
Mostra/amaga el codi
# Violin plot dels tres índexsdiv_long <- div_df %>%pivot_longer(cols =c(S, H_shannon, J_pielou),names_to ="index", values_to ="valor") %>%mutate(index =factor(index,levels =c("S", "H_shannon", "J_pielou"),labels =c("Riquesa (S)", "Shannon (H')", "Pielou (J')")))# Afegim etiquetes de significaciósig_df <-map_dfr(unique(div_long$index), function(idx) { sub <- div_long %>%filter(index == idx) amp <- sub %>%filter(proteccio =="AMP") %>%pull(valor) np <- sub %>%filter(proteccio =="NP") %>%pull(valor) p <-tryCatch(wilcox.test(amp, np, exact =FALSE)$p.value, error = \(e) NA)tibble(index = idx,y_pos =max(sub$valor, na.rm =TRUE) *1.05,etiq =sig_label(p))})p_violin <-ggplot(div_long, aes(x = proteccio, y = valor, fill = proteccio)) +geom_violin(alpha =0.6, trim =FALSE) +geom_boxplot(width =0.15, outlier.shape =NA, alpha =0.8) +geom_jitter(width =0.06, size =1.8, alpha =0.7) +geom_text(data = sig_df, aes(x =1.5, y = y_pos, label = etiq),size =5, inherit.aes =FALSE) +facet_wrap(~index, scales ="free_y") +scale_fill_manual(values = colors_prot) +labs(title ="Comparació dels índexs de diversitat (AMP vs NP)",subtitle ="* p<0.05 | ** p<0.01 | *** p<0.001 | ns = no significatiu",x =NULL, y ="Valor") +theme_bw(base_size =13) +theme(legend.position ="none")print(p_violin)
Mostra/amaga el codi
desa_grafic("02d_violin_diversitat.png", p_violin, ample =12, alt =5)
4 Corba d’acumulació d’espècies
Què és la corba d’acumulació?
La corba d’acumulació mostra com augmenta el nombre d’espècies a mesura que afegim més transsectes (mostrejos). Una corba que s’estabilitza indica que el mostreig és representatiu. S’usa la rarefacció (aleatorització de l’ordre dels transsectes, 999 vegades) per obtenir la corba esperada i la variabilitat (±SD).
Cada espècie és un punt. L’eix X ordena les espècies de més a menys abundants (rang). L’eix Y (escala logarítmica) mostra l’abundància relativa. Una corba molt pendent → poca equitativitat (unes poques espècies dominen). Una corba plana → comunitat equilibrada.
desa_grafic("04_rang_abundancia.png", p_rang, ample =10, alt =6)
Mostra/amaga el codi
rang_df %>%group_by(proteccio) %>%slice_min(rang, n =10) %>%select(Protecció = proteccio, Rang = rang, Espècie = especie, N_total, `Abund.rel.(%)`= abund_rel) %>%mutate(`Abund.rel.(%)`=round(`Abund.rel.(%)`, 2)) %>%kable(caption ="Top 10 espècies més abundants per grup de protecció") %>%kable_styling(bootstrap_options =c("striped", "hover"))
Top 10 espècies més abundants per grup de protecció
Protecció
Rang
Espècie
N_total
Abund.rel.(%)
AMP
1
Sarpa salpa
3840
39.49
AMP
2
Chromis chromis
1847
18.99
AMP
3
Diplodus vulgaris
982
10.10
AMP
4
Coris julis
867
8.92
AMP
5
Diplodus sargus
594
6.11
AMP
6
Symphodus tinca
588
6.05
AMP
7
Thalassoma pavo
216
2.22
AMP
8
Serranus cabrilla
214
2.20
AMP
9
Chelon labrosus
168
1.73
AMP
10
Serranus scriba
129
1.33
NP
1
Sarpa salpa
1449
40.68
NP
2
Chromis chromis
850
23.86
NP
3
Diplodus vulgaris
361
10.13
NP
4
Diplodus sargus
180
5.05
NP
5
Thalassoma pavo
135
3.79
NP
6
Coris julis
134
3.76
NP
7
Symphodus tinca
116
3.26
NP
8
Boops boops
101
2.84
NP
9
Chelon labrosus
100
2.81
NP
10
Serranus cabrilla
56
1.57
6 Abundància (densitat)
6.1 Totes les espècies: AMP i NP (gràfics separats)
Consell
Les barres mostren la densitat mitjana per espècie. Les barres d’error representen ±1 SE (error estàndard).
Mostra/amaga el codi
abund_trans <- dades %>%group_by(proteccio, estacio, transsecte, especie) %>%summarise(N =sum(nombre, na.rm =TRUE), .groups ="drop") %>%mutate(densitat = N /250)abund_resum <- abund_trans %>%group_by(proteccio, especie) %>%summarise(mitjana =mean(densitat, na.rm =TRUE),se =sd(densitat, na.rm =TRUE) /sqrt(n()),.groups ="drop" )ord_esp <- abund_resum %>%group_by(especie) %>%summarise(mg =mean(mitjana)) %>%arrange(desc(mg)) %>%pull(especie)abund_resum <- abund_resum %>%mutate(especie =factor(especie, levels = ord_esp))# Funció per fer el gràfic per un grupgrafic_abund <-function(prot, color_fill) { abund_resum %>%filter(proteccio == prot) %>%ggplot(aes(x = especie, y = mitjana)) +geom_col(fill = color_fill, alpha =0.85) +geom_errorbar(aes(ymin = mitjana - se, ymax = mitjana + se), width =0.4) +labs(title =paste("Densitat mitjana (±SE) per espècie –", prot),subtitle ="Tècnica: Comptatge visual",x =NULL, y ="Densitat (ind/250m²)") +theme_bw(base_size =11) +theme(axis.text.x =element_text(angle =45, hjust =1, face ="italic"),panel.grid.major.x =element_blank())}p_abund_AMP <-grafic_abund("AMP", "#FF7043")p_abund_NP <-grafic_abund("NP", "#2196F3")print(p_abund_AMP)
Mostra/amaga el codi
print(p_abund_NP)
Mostra/amaga el codi
desa_grafic("05a_densitat_global_AMP.png", p_abund_AMP, ample =12, alt =6)desa_grafic("05b_densitat_global_NP.png", p_abund_NP, ample =12, alt =6)
6.2 Per espècie individual (AMP vs NP en el mateix gràfic)
Consell
Selecciona les espècies que t’interessin canviant la llista especies_sel al bloc de codi. Per defecte, es mostra el conjunt d’espècies amb densitat > 0 en algun grup.
Mostra/amaga el codi
# Llistat d'espècies presentsespecies_all <-sort(unique(abund_trans$especie))# Gràfic per a cada espècieplots_abund_sp <-map(especies_all, function(sp) { df <- abund_trans %>%filter(especie == sp) kw <-tryCatch(kruskal.test(densitat ~ proteccio, data = df), error = \(e) NULL) sub <-if (!is.null(kw) &&!is.na(kw$p.value))paste0("KW: χ²=", round(kw$statistic, 2), ", p=", round(kw$p.value, 4)," ", sig_label(kw$p.value))else"Test no disponible" df %>%group_by(proteccio) %>%summarise(m =mean(densitat, na.rm =TRUE),se =sd(densitat, na.rm =TRUE) /sqrt(n()), .groups ="drop") %>%ggplot(aes(x = proteccio, y = m, fill = proteccio)) +geom_col(alpha =0.85, width =0.5) +geom_errorbar(aes(ymin = m - se, ymax = m + se), width =0.15, linewidth =0.8) +scale_fill_manual(values = colors_prot) +labs(title = sp, subtitle = sub, x =NULL, y ="Densitat (ind/250m²)") +theme_bw(base_size =11) +theme(legend.position ="none",plot.title =element_text(face ="italic", size =11))})names(plots_abund_sp) <- especies_all# Mostrem tots els gràfics i els desemfor (sp in especies_all) {print(plots_abund_sp[[sp]]) nom <-paste0("05c_densitat_", gsub(" ", "_", sp), ".png")desa_grafic(nom, plots_abund_sp[[sp]], ample =5, alt =4)}
6.3 Anàlisi Kruskal-Wallis per espècie
Nota
El test de Kruskal-Wallis és l’equivalent no paramètric de l’ANOVA. Compara les distribucions de densitat entre AMP i NP. Un p-valor < 0.05 indica diferències significatives.
P = Petits (talla < 1/3 de la talla màxima de l’espècie)
M = Mitjans (talla entre 1/3 i 2/3)
G = Grans (talla > 2/3)
La distribució de talles ens informa sobre l’estructura d’edats de la població i l’efecte de la pesca. En zones protegides esperem trobar més individus grans.
Mostra/amaga el codi
# Preparem les dades de tallestalles_base <- dades %>%filter(!is.na(talla)) %>%group_by(proteccio, transsecte, especie, talla) %>%summarise(N =sum(nombre, na.rm =TRUE), .groups ="drop")# Completem amb zeros explícits# SOLUCIÓ: complete() sense group_by previ, especificant totes les combinacionstalles_completa <- talles_base %>%ungroup() %>%complete(nesting(proteccio, transsecte, especie), talla,fill =list(N =0) )# Comprovaciócat("Files talles_base: ", nrow(talles_base), "\n")
PCA (Anàlisi de Components Principals): redueix la dimensionalitat de la matriu d’espècies (biomassa per estació). Estacions pròximes al biplot comparteixen composició similar. S’aplica la transformació de Hellinger abans del PCA.
RDA (Anàlisi de Redundància): combina la matriu d’espècies amb variables ambientals (profunditat, rugositat, cobertures). Permet identificar quins factors ambientals expliquen les diferències entre comunitats.
