Linking fish physiological biomarkers to habitat integrity in the Atlantic Forest

Biological Conservation
OR
Neotropical Ichthyology

Authors
Affiliations

Adamastor Coutinho Pinto

Programa de Pós-Graduação em Ecologia e Conservação (PPGEC), Universidade Estadual da Paraíba - Campus I, Rua Baraúnas, 351, Bairro Bodocongó, CEP 58109-753, Campina Grande, PB, Brazil.

Larissa Rafaela Caetano da Silva

Laboratório de Ecofisiologia Animal (LEFA), Departamento de Ciências Biológicas, Universidade Estadual da Paraíba – UEPB, Campus V, Rua Horácio Trajano de Oliveira, S/N, Cristo, CEP 58071-470, João Pessoa, PB, Brazil.

Mayara Mirelly da Silva Monteiro

Grupo de Ecologia de Rios do Semiárido, Laboratório de Ecologia, Departamento de Ciências Biológicas, Universidade Estadual da Paraíba – UEPB, Campus V, Rua Horácio Trajano de Oliveira, S/N, Cristo, CEP 58071-470, João Pessoa, PB, Brazil.

Alice da Silva Barros

Programa de Pós-Graduação em Ecologia e Conservação (PPGEC), Universidade Estadual da Paraíba - Campus I, Rua Baraúnas, 351, Bairro Bodocongó, CEP 58109-753, Campina Grande, PB, Brazil.

Enelise Marcelle Amado

Laboratório de Ecofisiologia Animal (LEFA), Departamento de Ciências Biológicas, Universidade Estadual da Paraíba – UEPB, Campus V, Rua Horácio Trajano de Oliveira, S/N, Cristo, CEP 58071-470, João Pessoa, PB, Brazil.

Elvio Sergio Figueredo Medeiros

Grupo de Ecologia de Rios do Semiárido, Laboratório de Ecologia, Departamento de Ciências Biológicas, Universidade Estadual da Paraíba – UEPB, Campus V, Rua Horácio Trajano de Oliveira, S/N, Cristo, CEP 58071-470, João Pessoa, PB, Brazil.

Published

August 30, 2026

Abstract
The Atlantic Forest is one of the world’s most important biodiversity hotspots. Despite occupying only a small fraction of its original distribution, it continues to support extraordinary levels of species richness and endemism. At the same time, it is one of the most threatened tropical ecosystems on Earth due to habitat loss, fragmentation, urban expansion, and agricultural development. The present study evaluated whether biomarkers of physiological stress in fish could detect environmental disturbances occurring within and around the Guaribas Biological Reserve (REBIO Guaribas), a strictly protected area of Atlantic Forest in northeastern Brazil. Despite sites within the administrative boundaries of the REBIO showing no signs alteration, from both environmental variables and biomarkers, significant differences in biomarker responses were observed among sampling sites, indicating that environmental conditions are not homogeneous throughout the stream studied. Elevated levels of all biomarkers at sites located outside or near the boundaries of the protected area suggest exposure to environmental stressors potentially associated with anthropogenic activities. These findings support the use of biochemical and physiological biomarkers as sensitive tools for detecting sublethal effects of environmental degradation, even in landscapes that contain legally protected areas.
Keywords

Ecological Health, Physiological Stress, REBIO Guaribas

1 Introduction

The Atlantic Forest is one of the world’s most important biodiversity hotspots (Myers et al., 2000). Despite occupying only a small fraction of its original distribution, it continues to support extraordinary levels of species richness and endemism (Marques et al., 2021). At the same time, it is one of the most threatened tropical ecosystems on Earth due to habitat loss, fragmentation, urban expansion, and agricultural development (Marques & Grelle, 2021). Conservation efforts in the Atlantic Forest have traditionally focused on terrestrial organisms such as plants, birds, and mammals. However, freshwater ecosystems are equally dependent on forest integrity and are often overlooked in conservation planning (Abilhoa et al., 2011; Gouveia et al., 2017). Furthermore, the predominance of small fish inhabiting small basins can lead to higher levels of endemism, resulting from allopatric speciation generated by low dispersion made by these fish and the consequent isolation of populations (Gouveia et al., 2017).

Padial et al. (2021) mention a wide range of threats to aquatic systems that can affect the Atlantic Forest, ranging from dam construction and flow alterations to the emerging threat of microplastic contamination (Eerkes-Medrano & Thompson, 2018). Pollution, from pesticides, fertilizers, heavy metals, pharmaceuticals or microplastics (Evans, 1987; Harmon, 2018), come mostly from urban and agricultural sources, and untreated urban sewage is amongst the most common (Costa et al., 2020). In agricultural areas, fish farming, ⁠agricultural runoff, and the ⁠destruction of riparian forests, remain as important threats to aquatic ecosystem integrity (Padial et al., 2021).

Ecological responses to pollution, often become visible only after environmental impacts have already affected freshwater organisms (Macêdo et al., 2019). Furthermore, the toxicity responsible for behavioral change depends on the organism, time of exposure and particle type, which demonstrates the need for physiological studies, specially in natural environments (Amoatey & Baawain, 2019). For these reasons, physiological biomarkers can provide particularly valuable early-warning indication (Ondei et al., 2020), as physiological, morphological and behavioural changes of species are affected in the ecosystem.

Fish are among the most widely used organisms to evaluate the health of aquatic systems, since their biochemical and physiological changes serve as biomarkers of environmental pollution (Cazenave et al., 2009; Ondei et al., 2020). One of the main physiological functions affected by pollutants is the osmoregulation, and its associated pH control. Both are predominantly done by gills, and are critical to gas exchange, ionic regulation and nitrogen excretion (Evans, 1987). The effects of this action on fish can be measured through a few blood markers, mainly lower levels of Chloride and Sodium ions, obtained by freshwater fish through exchange with internal ions (Evans, 1987), and higher levels of plasmatic proteins (Sabae & Mohamed, 2015). This osmoregulation impairment and consequent decrease in plasma ionic concentration can represent an osmotic challenge for tissues. The hydration level of gills and muscles is another physiological parameter associated with an animal’s capacity to cope with environmental pollution and is therefore subject to measurable changes induced by toxic compounds present in the environment. This active control was suggested by Ayrapetyan (2012) as a universal biomarker for detecting pollution in the environment (David et al., 2018; Sabae & Mohamed, 2015).

Another relevant biomarker is the abnormal activity of the multixenobiotic resistance (MXR) phenotype (Bard, 2000; Žaja et al., 2006), an innate gene of many aquatic species that is expressed in organs such as gills, kidneys, blood-brain barrier, liver and pancreas. It is activated by exposure to pollutants and results in the expression of P-glycoproteins (P-gp) in the membrane, which act to remove endogenous or exogenous toxic compounds out of the cell, protecting DNA and inducing toxic compounds excretion (Smital & Kurelec, 1998). The MXR activity represents the first line of biological defense, being induced by the presence of xenobiotics in the environment (Bard, 2000; Epel, 1998; Kurelec, 1992; Macêdo et al., 2019). Nevertheless, although the MXR activity explains how some species are more resistant to the presence of pollutants, some substances can reverse its effectiveness, being called chemosensitizers (or simply MXR-inhibitors) (Smital et al., 2004). Those compounds can interact with the P-gp activity, altering the mechanism regulation, and saturate or even reverse MXR activity. By saturating the P-gp pump, these chemosentizers allow other xenobiotics to enter the cell and increase its toxicity, even if the chemosensitizers themselves are not toxic. Precisely because they are mostly organic, natural, primarily harmless and very present in the materials emitted anthropogenic in nature, its detection should also be a priority in environmental risk studies (Kurelec et al., 2000).

In the present study, the hypothesis of markers linked to physiological responses of fish providing evidence of the environmental degradation was tested along a headwater stream that originates withing a protected area - Biological Reserve, or REBIO, according to BRASIL (2000) - and its course extends beyond its administrative boundaries (MMA, 2003). Given the importance of the Atlantic Forest to the biodiversity (Pereira et al., 2026), we aim with this study, (1) to identify physiological biomarkers in fish from a headwater stream, and (2) to associate these biomarkers to the ecological health of sites inside and surrounding a protected area of Atlantic Forest. This will enable decision makers to take action on the influence of anthropogenic threats and the ecological consequence of pollution exposure of the buffer zone.

2 Material and Methods

Study area and sampling design. This study was conducted in the Guaribas Biological Reserve (REBIO Guaribas), a protected area of Atlantic Forest in in northeastern Brazil. The reserve is located in the municipality of Mamanguape, state of Paraíba, in a heterogeneous landscape influenced by the adjacent Caatinga ecoregion. It was created in 1990 to preserve one of the last remnants of the Atlantic Forest in the region (MMA, 2003). Despite a REBIO being one of the most restricted types of protected areas (BRASIL, 2000), there are a number human settlements in its buffer zone. The management plan for the REBIO states several rural communities directly surrounding the reserve (Caiana, Pepina, Imbiribeira, Brejinho, Água Fria, João Pereira and Piabuçu), in addition to indigenous human populations (mostly the Potiguaras) (MMA, 2003). These populations are important because of their historic patterns of ocupation and land use, but also because their agricultural activities, fishing, modification of the natural vegetation for livestock and overal use of the natural resources from the buffer zone directly affect the conservation of the REBIO (Arruda et al., 2013; MMA, 2026). This has been aggravated in recent years by large scale use of water and land for cultivation of monocultures, like sugarcane (Heinrichs et al., 2017). The surrounding landscape is predominantly agricultural with other monocultures like coconut plantations , as well as the natural vegetation being modified for livestock. Other land uses include, farming, orchards and subsistence agriculture. These land uses can act as sources of diffuse pollution, especially through pesticide runoff, direct domestic effluents or contamination of underground water from septic tanks, creating a spatial gradient of environmental disturbance (MMA, 2003; Soler, 2004). Thus, the landscape is mainly rural and the land use causes direct impact to aquatic systems, by causing siltation of streams and contaminating groundwater through the use of septic tanks and releasing agricultural chemicals into the soil.

The climate in the region is classified as As’, characterized by a rainy season from February to July and a dry period from October to December (Peel et al., 2007). The average annual temperature varies from 24 to 26 °C, while the annual rainfall varies between 1,750 and 2,000 mm. The hydrographic network is composed of small headwater streams fed by rain, notably the Barro Branco and Caiana streams, both tributaries of the Camaratuba River (Gouveia et al., 2017).

Sampling was conducted along the Barro Branco stream between November 2023 and February 2024, comprising six sampling sites. The sampling sites were distributed longitudinally to capture environmental variability and gradients of human influence. Sites 1 to 3 were located within the protected area, where more preserved conditions were expected, while sites 4 to 6 were situated in the buffer zone adjacent to areas with greater human influence (Figure 1).

For the sake of comparison, we classified all sites into four categories defined a priori (Figure 2). Reference sites or Type I sites are sites located within the protected area and characterized by relatively preserved conditions, including intact riparian vegetation and no nearby crop fields (at least 1 km distant). Type I sites also fall into water-quality criteria, that are: low to moderate average temperatures (≤ 28 °C), relatively high average dissolved oxygen concentrations (≥ 6 mg/L), free-flowing water (average ≥ 0.025 m/s) and low average turbidity (≤ 15 NTU). Type II sites are located within the protected area but show some degree of environmental alteration, such as proximity to crop fields or deviations from one or more of the water-quality or hydrological characteristics used to define Type I sites. Type III sites occur outside the protected area but retain relatively good water-quality conditions comparable to those of Type I sites, although they may differ in other environmental characteristics associated with their location outside the protected area. Type IV sites occur outside the protected area and show evidence of environmental disturbance, including deterioration of water-quality conditions, not meeting at least one of the water-quality criteria, and/or alteration of habitat characteristics relative to the reference sites. Sites type III and IV were intentionally chosen not to meet the riparian forest and crop proximity criteria. Since they fall outside of the REBIO they are not expected to have intact riparian forest or to be at least 1 km distant from a crop field. These categories were latter confirmed by in locu environmental data measurements.

The ecophysiological evaluation of human impact on the forest, was mainly guided by the fish species inventory for the REBIO made by Gouveia et al. (2017). This study lists 18 species from five different orders: Characiformes (12 species), Perciformes (3 species), Synbranchiformes (1 species), Siluriformes (1 species) e Cyprinodontiformes (1 species). From those, two are exotic - Oreochromis niloticus, second most farmed fish species in the world and tolerant to osmotic variation, and Poecilia reticulata, largely used in mosquito control and resistant to stressful environments (Miranda, 2012). The predominance of small fish inhabiting small basins is thought to lead to low dispersion and consequent isolation of populations (Gouveia et al., 2017), associated with the fact that most species in the study area are phylogenetic close (with the predominance of Characidae) and are widely distributed, offers an appropriate design to access the entire local ichthyofauna using a few species, namely the genuses Hemigrammus and Astyanax.

Environmental data collection. In order to quantitatively confirm the classification of site category types (I, II, III and IV), environmental characteristics of each site were assessed based on four sets of variables: (a) site morphology, (b) water quality, (c) sediment composition, and (d) marginal habitat structure. Site morphology was characterized by measuring stream width (cm) and depth (cm) from three randomly selected transects. The slope of the margins was classified as 1 for a declivity of up to 20°, 2 for 20 to 30°, 3 for a declivity between 30 and 50° and 4 for declivity of the margins greater than 50°. Catchment-scale variables, such as elevation and stream length, were obtained using a handheld GPS and satellite imagery. Water quality was evaluated through physical and chemical variables measured with portable equipment, including temperature (°C) and dissolved oxygen concentration (mg/L). Water turbidity was determined in the laboratory using a turbidimeter, based on water samples collected at each sampling site. Water velocity (m/s) was estimated using the float method by Maitland (1990). Sediment composition and habitat physical structure were assessed following protocols adapted by Medeiros et al. (2008) from Pusey et al. (2004) and Mugodo et al. (2006). These variables were estimated from 9 to 12 one-meter quadrats distributed along the stream margins (terrestrial–aquatic interface). Within each quadrat, the percentage cover of sediment types (mud, sand, cobbles, small gravel, large gravel, rocks, and bedrock) and littoral and subaquatic habitat structures was visually estimated. Habitat structure variables included macrophytes, marginal grasses, leaf litter, attached algae, filamentous algae, overhanging vegetation, submerged vegetation, small woody debris, large woody debris, and root masses. All values of the environmental variables are henceforth presented as averages for sampling sites.

Three 250 mL replicate samples of water from each site were taken and refrigerated for concentration analysis (μg L⁻¹) of ammonia, nitrite, nitrate, orthophosphate (soluble reactive phosphorus) and total phosphorus. Nutrient concentrations were obtained based on Standard Methods for the Examination of Water (1998) and Mackereth et al. (1978).

Fish collection and maintenance. Fish were collected between October and November 2023 along the Barro Branco Stream at six previously established sampling points, located along the gradient of the stream, following a methodology adapted from Gouveia et al. (2017). Sampling was carried out during the day using a beach seine net (4 m long, 1.5 m high, and 5 mm mesh) and a hand net (50 cm wide and 5 mm mesh) operated with sweeping movements across a streach of 50 m (Medeiros et al., 2010). Fish sampling was performed randomly, therefore the occurrence of the species at each sampling site and the sample size varied according to natural densities of fish populations.

Fish were transferred to plastic bags containing water from the collection site and atmospheric oxygen, and transported to the Laboratory of Animal Ecophysiology (LEFA) at the State University of Paraíba. Fish were allowed to acclimated for a few hours to minimise stress from handling and transportation, and the ecophysiological experiments occurred on the same day of catch (Macêdo et al., 2019).

Experimental design. For analyses related to osmoregulatory function, blood samples were obtained from 78 individuals, including 38 fish collected within the reserve and 40 in external areas. The fish were anesthetized with eugenol, and blood aliquots were collected by cardiac puncture using pipettes with heparinized tips. The samples were kept on ice until laboratory processing and subsequently centrifuged for plasma separation. Due to the reduced volume of plasma obtained as a consequence of the small body size of the specimens, chloride ion analysis could not be performed. Branchial and muscle tissues from the same individuals were collected to determine tissue moisture content. Plasma and tissue samples were stored at −20 °C until laboratory analysis was performed. The remaining 111 individuals, comprising 54 specimens from locations within the reserve and 57 from external areas, were used for evaluation of the multixenobiotic resistance (MXR) phenotype, following the methodologies described by David et al. (2018); Macêdo et al. (2019); Santos et al. (2017).

Plasmatic protein dosage. For the the determination of plasmatic protein dosage, the plasma samples were defrosted and an aliquot of 2 µL of each was diluted with 8 µL of distilled water (⅕ proportion) to perform the total protein dosage, where the dye was also diluted in ⅕ proportion (BioRad® protein assay) (Bradford, 1976). The intensity of each sample coloration was measured with a microplate reader (SpectraMax i3) at 595 nm wave length, and compared to a standard curve (BSA protein). The results were transferred to a spreadsheet, where the average value in milligrams of protein dosage was calculated for each sampling site.

Moisture content. For the analysis of moisture content, parameter related to the osmoregulatory function, fish tissue samples (gill and muscle) were individually stored in 2 mL eppendorf tubes (previously weighted), initially weighted with the tissue moist, and, after 24 hours inside a 60°C oven, weighted with the dry tissue and the percentual content of tissue water was calculated (David et al., 2018).

MXR phenotype activity. The multixenobiotic resistance activity was analysed through the rhodamine B accumulation assay adapted from Smital & Kurelec (1998). Rhodamine B is a substrate of P-glycoprotein, the molecular basis of the MXR phenotype. Thus, after exposing the fish to this substrate, the analysis of the amount of rhodamine accumulated in the animal’s tissues (mainly the gills) reflects the activity of the MXR phenotype: the greater the accumulation, the lower the activity, and the lower the accumulation, the greater the activity. Therefore, specimens of fish from each sampling site (three inside the Conservative Unity and three outside, total of six sampling sites) were transferred to plastic aquariums containing dechlorinated water and rhodamine B at 2,5 µM concentration, staying in this condition for one hour. The aquarium remained under constant aeration and protected from direct light. After exposition, the fish were anesthetized with eugenol and sacrificed for gill sampling. The tissue samples were allocated in Eppendorfs tubes and weighted in an analytical scale. The weight of moist tissue was obtained by subtracting the weight of the previously weighted empty eppendorf tubes.The tissues were then homogenized with 500 µL of distilled water and the supernatant was transferred in triplicate of 100 µL to a 96-well microplate. The fluorescence intensity of the supernatant (corresponding to the intracellular rhodamine B fluorescence in the tissue) was measured in the microplate reader (SectraMax i3) at 544 nm excitation and 590 nm emission. The fluorescence value was then normalized by the moist tissue weight in milligrammes (Macêdo et al., 2019).

Statistical analysis. Environmental variables were assessed for multivariate collinearity, excluding redundant variables and combining those that conveyed similar information to reduce data redundancy, following the approach described by Legendre & Legendre (1998).

Code: Importing and organizing environment data
#MXR23 - ORGANIZANDO DADOS----

dev.off() #apaga os graficos, se houver algum
rm(list=ls(all=TRUE)) #limpa a memória
cat("\014") #limpa o console
getwd()
library(openxlsx)

##CARREGANDO MATRIZES BRUTAS----

habitat <- read.xlsx("data/rebio23-habitat.xlsx",
                     rowNames = T,
                     colNames = T,
                     sheet = "ambiente",
                     rows = 2:20) #startRow = 2

habitat[is.na(habitat)] <- 0
habitat

###REMOVENDO COLUNAS ZERADAS DE HABITAT----

sum <- colSums(habitat)
sum
zero_sum <- names(which(colSums(habitat) == 0))
zero_sum #nomes das colunas zeradas
m_part_cols <- habitat[(colSums(habitat) != 0)] #em != a exclamação inverte o sentido
zero_sum2 <- names(which(colSums(m_part_cols) == 0))
zero_sum2 #nomes das colunas zeradas
sum<-colSums(m_part_cols)
sum

t_grps <- read.xlsx("data/rebio23-habitat.xlsx",
                    rowNames = T,
                    colNames = T,
                    sheet = "grupos",
                    rows = 2:20) #startRow = 2
t_grps


###SEPARANDO NUTRIENTES----

nutri <- grep("^n\\.", colnames(m_part_cols), value = TRUE)
nutri
nutrientes <- m_part_cols[, nutri, drop = FALSE]

m_hab <- m_part_cols[, !colnames(m_part_cols) %in% nutri, drop = FALSE]
colnames(m_hab)

##SALVANDO MATRIZES FINAIS----

write.table(m_hab, "m_hab.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
write.table(nutrientes, "t_nutri.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
write.table(t_grps, "t_grps.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)
t_nutri <- read.csv("t_nutri.csv",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)
m_hab <- read.csv("m_hab.csv",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)
Code: Correlations of environment data
#MXR23----

##ORGANIZANDO DADOS----

dev.off()
rm(list=ls(all=TRUE))
cat("\014")

t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)
m_hab <- read.csv("m_hab.csv",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)

##CORRELOGRAMA E REMOÇÃO DE VARIÁVEIS REDUNDANTES OU DESNECESSÁRIAS----

library(psych)
colnames(m_hab)

png("fig-h.hab_pairs.png")
pairs.panels(m_hab[,15:24],
             method = "pearson", # correlation method
             scale = FALSE, lm = FALSE,
             hist.col = "#00AFBB", pch = 19,
             density = TRUE,  # show density plots
             ellipses = TRUE, # show correlation ellipses
             alpha = 0.5)
dev.off()

cor <- cor(m_hab)
cor

library(corrplot)
png("fig-hab_corrplot.png")
corrplot(cor, method = "circle")
dev.off()

##DELETANDO COLINEARES OU INDESEJADAS----

sink(file = "colineares.txt", append = F, split = T)
colnames(m_hab)
del_cols <- c("m.Depth_max_cm", "m.Depth_mar_cm")
#del_cols <- grep("^n\\.", colnames(m_hab), value = TRUE)
m_hab_part <- m_hab[, !(colnames(m_hab) %in% del_cols)]

##SOMANDO REDUNDANTES----

m_hab_part$h.Algae <- m_hab_part$h.Algae_f + m_hab_part$h.Algae_a
m_hab_part <- m_hab_part[, !(colnames(m_hab_part)
                           %in% c("h.Algae_f", "h.Algae_a"))]
m_hab_part$h.Debris <- m_hab_part$h.Debris_l + m_hab_part$h.Debris_s
m_hab_part <- m_hab_part[, !(colnames(m_hab_part)
                         %in% c("h.Debris_l", "h.Debris_s"))]

colnames(m_hab_part)
m_hab_part
sink()

write.table(m_hab_part, "m_hab_part.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
m_hab_part <- read.csv("m_hab_part.csv",
                       sep = ";", dec = ".",
                       row.names = 1,
                       header = TRUE,
                       na.strings = NA)
Code: Correlations tests of environment data
colnames(m_hab)
library(PerformanceAnalytics)
chart.Correlation(m_hab[, c(15:24)], histogram=TRUE)
# ***p<0.001, **p<0.01, *p<0.05, .p<0.10

library(Hmisc)
mat <- as.matrix(m_hab[, 15:24])
res <- rcorr(mat, type = "pearson")
cor_mat <- res$r
p_mat   <- res$P

stars <- ifelse(p_mat < 0.001, "***",
         ifelse(p_mat < 0.01,  "**",
         ifelse(p_mat < 0.05,  "*", "")))

cor_star <- ifelse(p_mat < 0.05,
                   paste0(round(cor_mat, 2), stars),
                   NA)

cor_star <- matrix(
  cor_star,
  nrow = nrow(cor_mat),
  dimnames = dimnames(cor_mat)
)

cor_star[upper.tri(cor_star, diag = TRUE)] <- NA
cor_star <- noquote(cor_star)
print(cor_star, na.print = "")

These variables were subsequently subjected to Principal Component Analysis (PCA) to evaluate multivariate correlations among sites. Prior to analysis, site morphology and water quality variables were square root transformed, whereas sediment composition and the marginal habitat structure (which were measured as percentages) were arcsine square-root transformed after relativization by column total (McCune & Grace, 2002). PCA was performed using the FactoMineR package in R (Lê et al., 2008). All variables were centered and scaled to unit of variance.

Code: Organizing Environment data PCA
#PCA fviz package----
#browseURL("https://www.sthda.com/english/articles/31-principal-component-methods-in-r-practical-guide/112-pca-principal-component-analysis-essentials/")

#ORGANIZANDO DADOS----

dev.off()
rm(list=ls(all=TRUE))
cat("\014")

t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)
m_hab_part <- read.csv("m_hab_part.csv",
                       sep = ";", dec = ".",
                       row.names = 1,
                       header = TRUE,
                       na.strings = NA)

colnames(m_hab_part)

###RELATIVIZAÇÕES E TRANSFORMAÇÕES----
library(vegan)
library(tidyverse)

# Selecionando variáveis que começam com "w" ou "m"
m_wm <- m_hab_part[, grep("^(w|m)", colnames(m_hab_part))]

# Selecionando variáveis que começam com "s" ou "h"
m_sh <- m_hab_part[, grep("^(s|h)", colnames(m_hab_part))]

# Removendo os dois primeiros caracteres dos nomes das colunas
colnames(m_wm) <- sub("^..", "", colnames(m_wm))
colnames(m_sh) <- sub("^..", "", colnames(m_sh))

# Transformações
m_wm <- sqrt(m_wm)
m_sh <- asin(sqrt(decostand(m_sh, method = "total", MARGIN = 2)))

# Juntando tudo em uma única tabela
hab_trns <- cbind(m_sh, m_wm)

# Salvado m_hab_trns
write.table(hab_trns, "m_hab_trns.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
m_hab_trns <- read.table("m_hab_trns.csv",
                         sep = ";", dec = ".",
                         row.names = 1,
                         header = TRUE,
                         na.strings = NA)

To test the hypothesis that physiological biomarkers of fish provide evidence of environmental degradation, all biomarker results were compared with those from Site 1, which served as the reference site for this study. Site 1 was located deep within the protected area, distant from anthropogenic influences, and met the criteria for a Type I site. After ascertaining the non-normality of the study biomarkers, comparisons were performed using one-sided Wilcoxon rank-sum tests to test whether biomarker values at each site were significantly greater than the mean observed at the reference site (Ross & Willson, 2017; Zar, 2010). Since no relevant different was observed across fish species, they were pooled withing sites. All statistical analyses were conducted in the R statistical environment using RStudio (R Core Team, 2017; R Studio Team, 2022).

3 Results

Environmental variables and nutrient concentration. Water velocity ranged from 0.022 to 0.293 m/s, with the highest values recorded at site 5. Width varied from 1.2 to 7.2 m, while depth ranged from 13.7 to 44.2 cm. Waters were, on average, well oxygenated (from 5.4 to 8.1 mg/L), and cool, with average temperatures from 25.2 to 27.8 °C. Conductivity was high, between 129.0 and 163.3 μS/cm, and pH, neutral to slightly acidic (from 5.5 to 7.0). Turbidity was generally low, although higher values were observed at sites 5 and 6, reaching 17.5 and 29.6 NTU, respectively. Elevation ranged from 41.5 to 112.8 m a.s.l., while the slope reached 50° for some sites. Mud and sand were the main substrate components across all sampled sites, with average cover ranging from 30 to 100% and from 0 to 70%, respectively. Habitat structure varied along the stream course. The reference site (site 1) was characterized by high contributions of overhanging vegetation (90%) and exclusively muddy substrates. Site 2 exhibited similar characteristics, whereas site 3 showed a more reduced habitat complexity. Sites 4–6 displayed greater variation in habitat composition. Overal, the habitat structures contributing most were overhanging vegetation, debris, and root masses, although these contributions varied among sites. Aquatic macrophytes and littoral grass were particularly abundant at site 4, contributing on average 35.8% and 39.2%, respectively, whereas root masses represented the dominant habitat structure at site 5 (76.7%). Sand cover increased downstream, reaching 70% at site 6, where overhanging vegetation was markedly reduced (Table 1).

Code: Environment data table
#TABELA DE HABITAT----

m_hab_part <- read.csv("m_hab_part.csv",
                       sep = ";", dec = ".",
                       row.names = 1,
                       header = TRUE,
                       na.strings = NA)
t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

library(dplyr)
library(tidyr)
m_trab <- m_hab_part %>%
  rename_with(~ gsub("_", ".", .))
m <- m_trab %>%
  group_by(PontoN = t_grps$PontoN) %>%
  summarise(across(where(is.numeric),
                   list(mean = mean, min = min, max = max)),
            .groups = 'drop') %>%
  pivot_longer(
    cols = -c(PontoN),
    names_to = c("Variable", ".value"),
    names_sep = "_"
  )

m <- as.data.frame(m)
m_wide <- m %>%
  mutate(stat_string = ifelse(Variable == c("m.Vel.m.s"),
                              paste0(round(mean, 3), " (", round(min, 3), "-", round(max, 3), ")"),
                              paste0(round(mean, 1), " (", round(min, 1), "-", round(max, 1), ")"))) %>%
  unite("Location", PontoN, sep = "_") %>%
  select(Variable, Location, stat_string) %>%
  pivot_wider(names_from = Location, values_from = stat_string)

m_wide
m_wide <- as.data.frame(m_wide)
m_wide <- m_wide[, c("Variable", "Ponto5", "Ponto6", "Ponto7", "Ponto8", "Ponto9", "Ponto10")] #sequência das colunas (ou linhas)
m_wide

# Salvado m_wide
write.table(m_wide, "m_wide_hab.txt",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
m_wide_hab <- read.table("m_wide_hab.txt",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)

# Exportando dados para Excel
library(openxlsx)
#write.xlsx(m_wide, file = "tabela de habitat.xlsx", rowNames = FALSE)
wb <- loadWorkbook("tabela de habitat.xlsx")
writeData(wb, sheet = "Sheet 1", x = m_wide)
saveWorkbook(wb, "tabela de habitat.xlsx", overwrite = TRUE)

# Tabela final ajustada no Excel
library(readxl)
hab_gttable <- read_excel(
  path  = "tabela de habitat.xlsx",
  sheet = "tabela_final",
  range = "A1:G26") #ajustar para o tamanho da tabela final 
hab_gttable
library(gt)
hab_gttable <- gt(hab_gttable)
hab_gttable <- sub_missing(hab_gttable,
                      columns = everything(), missing_text = "")
#hab_gttable <- fmt_number(hab_gttable, columns = "F.O.(%)", decimals = 1)
hab_gttable
gtsave(hab_gttable, "gt-hab_gttable_xlsx.html")
saveRDS(hab_gttable, "gt-hab_gttable_xlsx.rds")
Code: Environment data variable summary
# Escolher sumário de uma variavel
m
var <- "w.Temp.C"
m[m$Variable == var, "mean"] #cada valor de var
summary(m[m$Variable == var, "mean"]) #sumário dos valores de var

# Escolher sumário de um grupo de variáveis do df m
vars <- unique(grep("^h\\.", m$Variable, value = TRUE))
summaries <- list() #criam uma lista vazia para guardar os sumários
# Loop para cada variável do grupo e guarda em summaries
for (var in vars) {
  summaries[[var]] <- summary(m[m$Variable == var, "mean"])
}
#var is a temporary variable used in the for loop to iterate through
#each variable name that starts with "h."

summaries
summary_table <- do.call(rbind, lapply(summaries, as.data.frame.list))
round(summary_table, 2)
#sink(file = "summary_h.txt", split = TRUE)
round(summary_table[order(summary_table$Mean, decreasing = FALSE), ], 2)
#sink()

# Tabela limpa
#summary_table <- cbind(Variable = rownames(summary_table), summary_table)
#rownames(summary_table) <- NULL
#colnames(summary_table) <- c("Variable", "Min", "Q1", "Median", "Mean", "Q3", "Max")
#summary_table

Principal Component Analysis (PCA) described the overall environmental structure of the sampling sites and identified the main variables responsible for their spatial differentiation in terms of physicochemical characteristics, channel morphology, substrate composition, and marginal habitat structure (Figure 3). The first two principal components explained 57% of the total variance, with the first axis accounting for 31% and the second axis for 26%. The first axis represented a gradient separating sites I and II from sites type IV. Sites type I and II were associated with higher dissolved oxygen concentrations, higher altitude, muddy substrates, and overhanging vegetation, while sites type IV were characterized by greater channel width, average depth, littoral grass, sandy substrate, and higher turbidity. The upper portion of the second axis of the PCA was associated with algae, macrophytes, and channel width, while the lower portion was related to submerged vegetation, root masses, water velocity, turbidity, and sandy substrates. The site type III occupied an intermediate position along the first axis but was clearly distinguished along the second axis by its association with submerged vegetation, root masses, and higher water velocity. Overall, PCA indicated a marked environmental gradient between the sampling sites, with sites IV differing from the upstream sites mainly due to channel morphology and substrate characteristics, while the local habitat structure and aquatic vegetation explained the separation of site type III from the other sites.

Code: Environment data PCA - fviz
#PCA fviz package----
#browseURL("https://www.sthda.com/english/articles/31-principal-component-methods-in-r-practical-guide/112-pca-principal-component-analysis-essentials/")

#ORGANIZANDO DADOS----

dev.off()
rm(list=ls(all=TRUE))
cat("\014")

t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)
m_hab_trns <- read.csv("m_hab_trns.csv",
                       sep = ";", dec = ".",
                       row.names = 1,
                       header = TRUE,
                       na.strings = NA)

#PCA----

library("FactoMineR")
library("factoextra")

pca <- PCA(m_hab_trns, scale.unit = TRUE, ncp = 5, graph = TRUE)
print(pca)

eig.val <- get_eigenvalue(pca)
eig.val

#sink(file = "cor_matrix.txt", append = FALSE)
pca$var$cor
#sink()

fviz_eig(pca, addlabels = TRUE, ylim = c(0,30))

var <- get_pca_var(pca)
var

fviz_pca_var(pca, col.var = "black")

library("corrplot")
corrplot(var$cos2, is.corr=FALSE)

fviz_cos2(pca, choice = "var", axes = 1:2)

fviz_pca_var(pca, col.var = "cos2",
             gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
             repel = TRUE # Avoid text overlapping
)

fviz_pca_var(pca, alpha.var = "cos2")

corrplot(var$contrib, is.corr=FALSE)

fviz_contrib(pca, choice = "var", axes = 1, top = 10)
fviz_contrib(pca, choice = "var", axes = 2, top = 10)
fviz_contrib(pca, choice = "var", axes = 1:2, top = 10)

fviz_pca_var(pca, col.var = "contrib",
             gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07")
)
fviz_pca_var(pca, alpha.var = "contrib")

##GROUPING BY KMEANS----

set.seed(123)
res.km <- kmeans(var$coord, centers = 3, nstart = 25)
grp <- as.factor(res.km$cluster)
# Color variables by groups
fviz_pca_var(pca, col.var = grp,
             palette = c("#0073C2FF", "#EFC000FF", "#868686FF"),
             legend.title = "Cluster")

ind <- get_pca_ind(pca)
ind
ind$contrib

##BIPLOTS----

fviz_pca_ind(pca)

fviz_pca_ind(pca, col.ind = "cos2",
             gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
             repel = TRUE #avoid text overlapping (slow if many points)
)

fviz_pca_ind(pca, pointsize = "cos2",
             pointshape = 21, fill = "#E7B800",
             repel = TRUE #avoid text overlapping (slow if many points)
)

fviz_pca_ind(pca, col.ind = "cos2", pointsize = "cos2",
             gradient.cols = c("#00AFBB", "#E7B800", "#FC4E07"),
             repel = TRUE #avoid text overlapping (slow if many points)
)

fviz_cos2(pca, choice = "ind")

fviz_contrib(pca, choice = "ind", axes = 1:2)

# Create a random continuous variable of length 23,
# Same length as the number of active individuals in the PCA
t_grps <- t_grps[rownames(pca$ind$coord), ]
set.seed(123)
my.cont.var <- t_grps$Ref_site
my.cont.var
# Color individuals by the continuous variable
fviz_pca_ind(
  pca,
  habillage = my.cont.var,
  palette = c("blue", "yellow", "red", "pink"),
  repel = FALSE
)

fviz_pca_ind(pca,
             geom.ind = "point", #show points only (nbut not "text")
             col.ind = as.factor(my.cont.var), #color by groups
             palette = "grey",
             addEllipses = TRUE, #concentration ellipses
             legend.title = "Groups"
)

length(unique(my.cont.var))
fviz_pca_ind(pca,
             geom.ind = "point",
             col.ind = t_grps$Ref_site,
             palette = "grey",
             addEllipses = TRUE, ellipse.type = "confidence",
             mean.point = TRUE,
             legend.title = "Groups")

fviz_pca_biplot(pca, repel = TRUE,
                col.var = "#2E9FDF", #variables color
                col.ind = "#696969"  #individuals color
)

fviz_pca_biplot(pca,
                col.ind = t_grps$Ref_site, palette = "jco",
                addEllipses = TRUE, elipse.type = "euclid",
                label = "var",
                col.var = "black", repel = TRUE,
                legend.title = "Species")

### FILTRANDO AS VARIÁVEIS
var_coord <- pca$var$coord[,1:2]
keep_vars <- apply(abs(var_coord) > 0.5,1,any)
vars_to_keep <- rownames(var_coord)[keep_vars]
pca_filt <- pca
pca_filt$var$coord    <- pca$var$coord[vars_to_keep,,drop=FALSE]
pca_filt$var$cos2     <- pca$var$cos2[vars_to_keep,,drop=FALSE]
pca_filt$var$contrib  <- pca$var$contrib[vars_to_keep,,drop=FALSE]

### AJUSTE MANUAL DAS COORDENADAS DOS VETORES
var_coords <- as.data.frame(pca_filt$var$coord)
var_coords$Dim.1 <- var_coords$Dim.1*4.2
var_coords$Dim.2 <- var_coords$Dim.2*4.2

## Caso queira editar manualmente
#fix(var_coords) #nomes das colunas
#write.table(var_coords,"coords_var.txt")

var_coords <- read.table("coords_var.txt",
                         header=TRUE,
                         row.names=1)

### AJUSTE MANUAL DOS INDIVÍDUOS
ind_coords <- as.data.frame(pca_filt$ind$coord)

## Caso queira editar manualmente
#fix(ind_coords) #nomes das linhas
#write.table(ind_coords,"coords_ind.txt")

ind_coords <- read.table("coords_ind.txt",
                         header=TRUE,
                         row.names=1)

############################################################
### GRÁFICO FINAL
############################################################

library(factoextra)
library(ggplot2)

t_grps$Ref_site <- as.factor(t_grps$Ref_site)

#png("fig-hab_PCA_fviz.png", width=11, height=9, units="in", res=400)

library(ggplot2)

ggplot() +
  
  geom_hline(yintercept = 0, linetype = 2) +
  geom_vline(xintercept = 0, linetype = 2) +
  
  geom_point(
    data = ind_coords,
    aes(Dim.1, Dim.2, shape = t_grps$Ref_site),
    size = 3,
    colour = "black"
  ) +
  
  geom_segment(
    data = var_coords,
    aes(x = 0, y = 0,
        xend = Dim.1,
        yend = Dim.2),
    colour = "red",
    arrow = arrow(length = unit(0.25,"cm"))
  ) +
  
  geom_text(
    data = transform(
      var_coords,
      nudge_y = ifelse(Dim.2 >= 0, 0.25, -0.25)
    ),
    aes(
      x = Dim.1,
      y = Dim.2 + nudge_y,
      label = rownames(var_coords)
    ),
    colour = "red",
    size = 4,
    nudge_x = 0.19
  ) +
  
  scale_shape_manual(values = c(15,16,0,1,2,17)) +
  
  labs(
    x = paste0("Dim1 (", round(pca$eig[1,2],1), "%)"),
    y = paste0("Dim2 (", round(pca$eig[2,2],1), "%)"),
    shape = "Sites"
  ) +
  
  coord_equal() +
  theme_minimal()

#dev.off()

Averaged nutrient concentrations varied among sampling sites (Table 2). Ammonia concentrations ranged from 54.8 µg L⁻¹ at the reference site (Site 1) to 243.7 µg L⁻¹ at Site 6, the highest value recorded during the study. Elevated ammonia concentrations were also observed at Sites 4 (81.6 µg L⁻¹) and 5 (80.8 µg L⁻¹), whereas Sites 1–3 exhibited comparatively lower concentrations (54.8–72.0 µg L⁻¹). Nitrite concentrations were too low to be detectable for at all sampling sites. Nitrate concentrations ranged from 5.5 to 18.8 µg L⁻¹, with the highest values recorded at Sites 4 and 5. Orthophosphate was detected only at Site 3 (1.0 µg L⁻¹), whereas total phosphorus ranged from undetectable levels at Site 1 to 50.1 µg L⁻¹ at Site 3, with relatively high concentrations also recorded at Site 6 (38.7 µg L⁻¹).

Code: Nutrient data table
#TABELA DE NUTRIENTES----

##ORGANIZANDO DADOS----

dev.off()
rm(list=ls(all=TRUE))
cat("\014")

t_nutri <- read.csv("t_nutri.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

t_grps <- read.csv("t_grps.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

nutri <- t_nutri

library(dplyr)
library(tidyr)
t_grps <- t_grps[seq_len(nrow(nutri)), , drop = FALSE]
m_trab <- nutri %>%
  rename_with(~ gsub("_", ".", .))
m <- m_trab %>%
  group_by(PontoN = t_grps$PontoN) %>%
  summarise(across(where(is.numeric),
                   list(mean = mean, min = min, max = max)),
            .groups = 'drop') %>%
  pivot_longer(
    cols = -c(PontoN),
    names_to = c("Variable", ".value"),
    names_sep = "_"
  )

m <- as.data.frame(m)
colnames(m_trab)
m_wide <- m %>%
  mutate(stat_string = ifelse(Variable == c("Ammonia.NH3.ug.L"),
                              paste0(round(mean, 1), " (", round(min, 1), "-", round(max, 1), ")"),
                              paste0(round(mean, 1), " (", round(min, 1), "-", round(max, 1), ")"))) %>%
  unite("Location", PontoN, sep = "_") %>%
  select(Variable, Location, stat_string) %>%
  pivot_wider(names_from = Location, values_from = stat_string)

m_wide
m_wide <- as.data.frame(m_wide)
m_wide <- m_wide[, c("Variable", "Ponto5", "Ponto6", "Ponto7", "Ponto8", "Ponto9", "Ponto10")] #sequência das colunas (ou linhas)
m_wide

# Salvado m_wide
write.table(m_wide, "m_wide_nutri.txt",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
m_wide_nutri <- read.table("m_wide_nutri.txt",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)

# Exportando dados para Excel
library(openxlsx)
#write.xlsx(m_wide, file = "tabela de nutrientes.xlsx", rowNames = FALSE)
wb <- loadWorkbook("tabela de nutrientes.xlsx")
writeData(wb, sheet = "Sheet 1", x = m_wide_nutri)
saveWorkbook(wb, "tabela de nutrientes.xlsx", overwrite = TRUE)

# Tabela final ajustada no Excel
library(readxl)
nutri_gttable <- read_excel(
  path  = "tabela de nutrientes.xlsx",
  sheet = "tabela_final",
  range = "A1:G6") #ajustar para o tamanho da tabela final 
nutri_gttable
library(gt)
nutri_gttable <- gt(nutri_gttable)
nutri_gttable <- sub_missing(nutri_gttable,
                      columns = everything(), missing_text = "")
nutri_gttable
gtsave(nutri_gttable, "gt-nutri_gttable_xlsx.html")
saveRDS(nutri_gttable, "gt-nutri_gttable_xlsx.rds")

Fish community. A total of six fish species were collected under a standardized sampling effort, so species representation reflected their natural abundances in the study area. Among the families recorded, Characidae was the richest, comprising five species (Astyanax bimaculatus, Astyanax fasciatus, Hemigrammus marginatus, Hemigrammus rodwayi, and Hemigrammus unilineatus), whereas Crenuchidae was represented by a single species, Characidium bimaculatum, from which only one specimen was collected. Hemigrammus unilineatus was the most abundant species, representing 66.3% of all individuals analyzed, followed by Astyanax fasciatus (14.8%) and Hemigrammus rodwayi (11.2%). Astyanax bimaculatus accounted for 4.6% of the specimens, whereas Hemigrammus marginatus and Characidium bimaculatum represented only 2.5% and 0.5%, respectively. Since species representation in the physiological analyses reflected their natural abundances, the number of individuals available was divided for each subsequent biomarker analysis, in order to reflect the natural composition of the fish assemblage. Of the total fish collected, 118 specimens were submitted to the rhodamine B accumulation assay and 78 were used for total plasma protein and moisture content analyses (Figure 4). Given its low count, Characidium bimaculatum was removed from all biomarker analyses.

Code: Counts of fish
dev.off() #apaga os graficos, se houver algum
rm(list=ls(all=TRUE)) #limpa a memória
cat("\014") #limpa o console

getwd()
library(openxlsx)
mxr23 <- read.xlsx("D:/Elvio/OneDrive/MSS/_rebio23-mxr/mxr23_Q/data/rebio23-peixes_5-10.xlsx",
                   rowNames = T,
                   colNames = T,
                   sheet = "peixes-fisio")
mxr23#[1:5,1:5] mostra apenas as linhas e colunas de 1 a 5.

###### ARRUMANDO OS DADOS #####

mxr23_part <- mxr23[!row.names(mxr23) %in% c("4-9-85", #NA
                                       "1-10-7"),] #As.bim grande

mxr23_part[c("1-10-4","1-10-5"), c("Rodamina_rfu_mg", "TH_musculo_p",
                             "TH_branquia_p", "Proteinas_Totais_mg"
                             )] <- NA
mxr23_part <- mxr23_part[mxr23_part$Especie != "Characidium bimaculatum", ]

write.table(mxr23_part, "mxr23_part.csv",
            sep = ";", dec = ".", #"\t",
            row.names = TRUE,
            quote = TRUE,
            append = FALSE)
data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

##### RESUMO DOS DADOS #####

library(dplyr)
library(tidyr)

grouped_data <- group_by(data, Especie, Amostra_tipo)
summarised_data <- summarise(grouped_data, count = n())
species_count_table <- spread(summarised_data, Amostra_tipo, count, fill = 0)
print(species_count_table)

species_count_table <- species_count_table %>%
  rowwise() %>%
  mutate(N = sum(c_across(where(is.numeric))))

resumo <- as.data.frame(species_count_table)
resumo

rownames(resumo) <- resumo[,1] #tem  que ser um df
resumo[,1] <- NULL

soma <- apply(resumo,2,sum)
soma

soma_row <- as.data.frame(t(soma))
rownames(soma_row) <- "Soma"

resumo <- rbind(resumo, soma_row)
resumo

###### RESUMO SPP POR UA #####

resumo2 <- data %>%
  group_by(Site = Site, Especie = Especie, Amostra_tipo = Amostra_tipo) %>%
  summarise(Count = n()) %>%
  arrange(Site, Especie, Amostra_tipo)

resumo2 <- as.data.frame(resumo2)
resumo2

write.table(resumo2, file = "resumo.csv", row.names = T, sep = "\t")
read.table("resumo.csv", check.names = F)

resumo2_wide <- resumo2 %>%
  pivot_wider(names_from = c(Site, Amostra_tipo), values_from = Count, values_fill = list(Count = 0))
resumo2_wide <- as.data.frame(resumo2_wide)
resumo2_wide

# Create the gt table
library(gt)
gt_table <- resumo2_wide %>%
  gt() %>%
  tab_header(
    title = "Species Count by Site and Sample Type"
  ) %>%
  cols_label(
    Especie = "Species"
  ) %>%
  fmt_number(
    columns = everything(),
    decimals = 0
  ) %>%
  cols_width(
    everything() ~ px(100)
  ) %>%
  opt_table_outline()

gt_table

##### GRÁFICO POR N #####

# Filtrar por Amostra_tipo "xx"
filtered_Amostra_tipo <- data %>%
  filter(Amostra_tipo != "xx")
data <- filtered_Amostra_tipo

# Count occurrences of each species by Amostra_tipo and Site
species_count <- data %>%
  group_by(Site, Especie, Amostra_tipo) %>%
  summarise(Count = n(), .groups = 'drop')

# Plot the counts
library(ggplot2)

counts <- ggplot(species_count, aes(x = Amostra_tipo, y = Count, fill = Especie)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  facet_wrap(
  ~ Site,
  labeller = labeller(
    Site = c(
      "B-P05" = "Site 1",
      "B-P06" = "Site 2",
      "B-P07" = "Site 3",
      "B-P08" = "Site 4",
      "B-P09" = "Site 5",
      "B-P10" = "Site 6"
    )
  )
) +
   scale_x_discrete(
    labels = c(
      "rd" = "MXR",
      "th" = "Moisture"
    )
  ) +
  scale_fill_grey(start = 0.2, end = 0.8) +
  labs(title = "",
       x = "Sample use",
       y = "Count",
       fill = "Species") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))
counts
ggsave(counts, dpi = 300, filename = "fig-counts.png", bg = "white")


# Count occurrences of each species by Amostra_tipo and Site
species_count <- data %>%
  group_by(Site, Especie, Amostra_tipo) %>%
  summarise(Count = n(), .groups = 'drop')
# Plot the counts with species on the x-axis
ggplot(species_count, aes(x = Especie, y = Count, fill = Amostra_tipo)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  facet_wrap(~ Site) +
  labs(title = "Count of Each Species per Amostra_tipo for Each Site",
       x = "Species",
       y = "Count") +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

Plasmatic protein dosage. The number of plasmatic protein dosages (n=44) was smaller than the total number of fish collected for this part of the analysis (n=78), because of the difficulty on collecting blood from smaller specimens and the need to collect a minimum amount of plasma (Figure 5).

Code: Total Plasmatic Proteins (mg)
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$Proteinas_Totais_mg))

##### GRAFICO POR ANÁLISE #####

#N #data <- data %>% mutate(N = 1)
#RPC
#Rodamina_rfu_mg
#TH_musculo_p
#TH_branquias_p
#Proteinas_Totais_mg

# Convert the column to numeric, replace commas with dots for proper conversion
data$Proteinas_Totais_mg <- as.numeric(gsub(",", ".", data$Proteinas_Totais_mg))
VAR <- deparse(substitute(data$Proteinas_Totais_mg))
var <- "Proteinas_Totais_mg" %>% rlang::sym()
str(data)

# Calculate the average for each combination of species, Amostra_tipo, and Site
average <- summarise(
  group_by(data, Site, Especie, Amostra_tipo),
  Average = mean(!!var, na.rm = TRUE),
  n = sum(!is.na(!!var)),
  SE = ifelse(
    sum(!is.na(!!var)) > 1,
    sd(!!var, na.rm = TRUE) / sqrt(sum(!is.na(!!var))),
    0
  ),
  .groups = "drop"
)
average <- average[
  average$Amostra_tipo == "th" & !is.nan(average$Average),
]
# Plot the average
ggplot(average, aes(x = Especie, y = Average, fill = Amostra_tipo)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  facet_wrap(~ Site) +
  labs(title = paste("Average POR ANÁLISE para", VAR, "per Species and Amostra_tipo for Each Site"),
       x = "Species",
       y = paste("Average", VAR)) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

##### PLOT FINAL #####

prot_tots <- 
  ggplot(average, aes(x = Especie, y = Average, fill = Amostra_tipo)) +
  geom_bar(stat = "identity", position = position_dodge(), fill = "grey70") +
    geom_errorbar(
    aes(
      ymin = Average - SE,
      ymax = Average + SE
    ),
    width = 0.2
  ) +
  geom_text(
    aes(
      y = 10,
      label = paste0("n=", n)
    ),
    size = 3
  ) +
  facet_wrap(
  ~ Site,
  labeller = labeller(
    Site = c(
      "B-P05" = "Site 1",
      "B-P06" = "Site 2",
      "B-P07" = "Site 3",
      "B-P08" = "Site 4",
      "B-P09" = "Site 5",
      "B-P10" = "Site 6"
    )
  )
) +
  scale_fill_grey(start = 0.2, end = 0.8) +
  labs(title = "",
       x = "Species",
       y = "Plasmatic protein dosage (mg)",
       fill = "Species") +
  theme_grey() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        legend.position = "none")
prot_tots
ggsave(prot_tots, dpi = 300, filename = "fig-prot_tots2.png", bg = "white")
Code: Total Plasmatic Proteins (mg) per species
##### GRÁFICO POR ESPÉCIE #####

# Filter out the species Astyanax bimaculatus and Astyanax fasciatus

filtered_data <- data

filtered_data <- data %>%
  filter(!Especie %in% c("Astyanax bimaculatus",
                        "Astyanax fasciatus",
                        "Characidium bimaculatum")) #EXCLUI USANDO O "!"
filtered_data <- data %>%
  filter(Especie %in% c("Hemigrammus unilineatus")) #ESCOLHE

# Calculate the average for each combination of Amostra_tipo, Especie, and Site
average <- filtered_data %>%
  group_by(Site, Amostra_tipo, Especie) %>%
  summarise(Average = mean(!!var, na.rm = TRUE), .groups = 'drop')
# Plot the average RPC with Amostra_tipo on the x-axis
ggplot(average, aes(x = Amostra_tipo, y = Average, fill = Especie)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  facet_wrap(~ Site) +
  labs(title = paste("Average", VAR, "per Amostra_tipo for Each Site"),
       x = "Amostra Tipo",
       y = paste("Average", VAR)) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

##### CORRELAÇÕES #####

data <- filtered_data

colnames(data)
cor <- data %>%
  select(Ponto, Ref_site, Compr.total_mm, Compr.padrao_mm, Altura_mm,
  Peso_g, RPC, Rodamina_rfu_mg, TH_branquia_p, TH_musculo_p,
  Proteinas_Totais_mg)
str(cor)
cor[] <- lapply(cor, function(x) as.numeric(gsub(",", ".", x)))

plot(cor)
plot(cor[,3:7])

cor(cor,method="pearson")
cor(cor[,1:3], method="spearman")

library(psych)
pairs.panels(cor[,8:11],
             method = "pearson", # correlation method
             scale = FALSE, lm = FALSE,
             hist.col = "#00AFBB", pch = 19,
             density = TRUE,  # show density plots
             ellipses = TRUE, # show correlation ellipses
             alpha = 0.5
)

cor %>%
  with(cor.test(RPC, Peso_g, method = "pearson"))

with(cor, plot(scale(RPC), scale(Peso_g)))
with(cor, plot(RPC, Peso_g))

plot((m_part$"Overh_veg" - mean(m_part$"Overh_veg")) / sd(m_part$"Overh_veg"))
#plot((m_trns$"m.elev" - mean(m_trns$"m.elev")) / sd(m_trns$"m.elev"))
with(cor, plot((RPC - mean(Peso_g)) / sd(Peso_g)))
with(cor, plot((RPC - mean(RPC)) / sd(RPC),
                  (Peso_g - mean(Peso_g)) / sd(Peso_g)))

Plasma total protein concentrations differed between sampling sites when compared with the reference site (Site 1). Mean protein concentrations were significantly higher at the sites 4, 5 and 6 (91.25 ± 5.68, 74.22 ± 3.31 and 74.79 ± 6.55 mg, respectively) than the reference site (33.53 ± 0.32 mg; Wilcoxon rank-sum test, p < 0.02). In contrast, no significant increase was observed at Sites 2 and 3 (29.41 ± 1.59 and 33.75 ± 0.80 mg, respectively; p > 0.4). The highest mean protein concentration was recorded at Site 4, representing an approximately 2.7-fold increase relative to the reference site.

Code: Total Plasmatic Proteins t-tests: one-sided two-sample t-test
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$Proteinas_Totais_mg))
var <- "Proteinas_Totais_mg"

# Reference site
ref <- data[[var]][
  data$Site == "B-P05" &
  !is.na(data[[var]])
]

# Sites to compare
sites <- setdiff(unique(data$Site), "B-P05")

# Run tests
results <- data.frame()

for(site in sites){

  test_data <- data[[var]][
  data$Site == site &
  !is.na(data[[var]])
]

  tt <- t.test(
    test_data,
    ref,
    alternative = "greater"
  )

  results <- rbind(
  results,
  data.frame(
    Site = site,
    Mean_Site = mean(test_data),
    SD_Site = sd(test_data),
    SE_Site = sd(test_data) / sqrt(length(test_data)),
    Mean_Ref = mean(ref),
    SD_Ref = sd(ref),
    SE_Ref = sd(ref) / sqrt(length(ref)),
    N_Site = length(test_data),
    N_Ref = length(ref),
    t = unname(tt$statistic),
    df = unname(tt$parameter),
    p = tt$p.value
  )
)
}

results

results$Significant <- ifelse(results$p < 0.05, "Yes", "No")

results$p <- ifelse(
  results$p < 0.001,
  "<0.001",
  sprintf("%.3f", results$p)
)
#sink("t-tests_prot_tots.txt", append = TRUE)
results
#sink()

##### NÃO-PARAMETRICO #####

results_w <- data.frame()

for(site in sites){

  test_data <- data[[var]][
    data$Site == site &
    !is.na(data[[var]])
  ]

  wt <- wilcox.test(
    test_data,
    ref,
    alternative = "greater"
  )

  results_w <- rbind(
    results_w,
    data.frame(
      Site = site,
      W = unname(wt$statistic),
      p = wt$p.value
    )
  )
}

results_w$Significant <- ifelse(results_w$p < 0.05, "Yes", "No")

results_w$p <- ifelse(
  results_w$p < 0.001,
  "<0.001",
  sprintf("%.3f", results_w$p)
)

#sink("t-tests_np.txt", append = FALSE)
results_w
#sink()

Moisture content. Gill and muscle tissues moisture content showed similar results to plasmatic protein dosage, with gill moisture concentrations being significantly higher than those observed at the reference site (Site 1; 82.82 ± 1.75%) for sites 5 and 6 (mean values of 90.13 ± 1.19% and 93.89 ± 0.98%, Wilcoxon rank-sum test; p < 0.004). No significant increases were detected at Sites 2, 3, or 4 (84.35 ± 1.60%, 85.63 ± 1.94%, and 82.84 ± 0.75%, respectively; p > 0.1 in all cases).

The muscle moisture concentrations at the reference site averaged 80.31 ± 1.48%. Significantly higher values were observed at Sites 4, 5 and 6, where mean moisture concentrations reached 82.38 ± 0.62%, 84.29 ± 0.88% and 86.02 ± 2.49%, respectively (Wilcoxon rank-sum test; p ≤ 0.04). In contrast, moisture concentrations at Sites 2 and 3 (80.11 ± 0.66% and 80.58 ± 0.94%) did not differ significantly from those at the reference site (p > 0.39) (Figure 6).

Code: Moisture content (gills and muscle)
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$TH_branquia_p))
table(data$Especie, is.na(data$TH_musculo_p))

##### GRAFICO POR ANÁLISE #####

#N #data <- data %>% mutate(N = 1)
#RPC
#Rodamina_rfu_mg
#TH_musculo_p
#TH_branquias_p
#Proteinas_Totais_mg

library(dplyr)
library(tidyr)
library(ggplot2)

# Convert the column to numeric, replace commas with dots for proper conversion
data$TH_branquia_p <- as.numeric(gsub(",", ".", data$TH_branquia_p))
data$TH_musculo_p  <- as.numeric(gsub(",", ".", data$TH_musculo_p))
str(data)

# Long format
th_data <- pivot_longer(
  data,
  cols = c(TH_branquia_p, TH_musculo_p),
  names_to = "Tecido",
  values_to = "TH"
)

# Calculate the average for each combination of species, Amostra_tipo, and Site
average <- summarise(
  group_by(th_data, Site, Especie, Tecido),
  Average = mean(TH, na.rm = TRUE),
  n = sum(!is.na(TH)),
  SE = ifelse(
    n > 1,
    sd(TH, na.rm = TRUE) / sqrt(n),
    0
  ),
  .groups = "drop"
)
# Remove groups with no observations
average <- average[average$n > 0, ]

##### PLOT FINAL #####

prot_th <-
  ggplot(
    average,
    aes(
      x = Especie,
      y = Average,
      fill = Tecido
    )
  ) +
  geom_col(
    position = position_dodge(width = 0.9)
  ) +
  geom_errorbar(
    aes(
      ymin = Average - SE,
      ymax = Average + SE
    ),
    width = 0.2,
    position = position_dodge(width = 0.9)
  ) +
  geom_text(
    aes(
      y = 5,
      label = paste0("", n)
    ),
    size = 3,
    position = position_dodge(width = 0.9)
  ) +
  facet_wrap(
    ~Site,
    labeller = labeller(
      Site = c(
        "B-P05" = "Site 1",
        "B-P06" = "Site 2",
        "B-P07" = "Site 3",
        "B-P08" = "Site 4",
        "B-P09" = "Site 5",
        "B-P10" = "Site 6"
      )
    )
  ) +
  scale_fill_grey(
    start = 0.4,
    end = 0.8,
    labels = c(
      "TH_branquia_p" = "Gill",
      "TH_musculo_p" = "Muscle"
    )
  ) +
  labs(
    x = "Species",
    y = "Moisture concentration (%)",
    fill = "Tissue"
  ) +
  theme_grey() +
  theme(
    axis.text.x = element_text(
      angle = 45,
      hjust = 1
    ),
    strip.background = element_rect(
      fill = "grey85",
      colour = "grey85"
    )
  )
prot_th
ggsave(prot_th, dpi = 300, filename = "fig-prot_th.png", bg = "white")
Code: Moisture content (gills and muscle) t-tests: one-sided two-sample t-test
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$TH_branquia_p))
table(data$Especie, is.na(data$TH_musculo_p))

var <- "TH_musculo_p"

# Reference site
ref <- data[[var]][
  data$Site == "B-P05" &
  !is.na(data[[var]])
]

# Sites to compare
sites <- setdiff(unique(data$Site), "B-P05")

# Run tests
results <- data.frame()

for(site in sites){

  test_data <- data[[var]][
  data$Site == site &
  !is.na(data[[var]])
]

  tt <- t.test(
    test_data,
    ref,
    alternative = "greater"
  )

  results <- rbind(
  results,
  data.frame(
    Site = site,
    Mean_Site = mean(test_data),
    SD_Site = sd(test_data),
    SE_Site = sd(test_data) / sqrt(length(test_data)),
    Mean_Ref = mean(ref),
    SD_Ref = sd(ref),
    SE_Ref = sd(ref) / sqrt(length(ref)),
    N_Site = length(test_data),
    N_Ref = length(ref),
    t = unname(tt$statistic),
    df = unname(tt$parameter),
    p = tt$p.value
  )
)
}

results

results$Significant <- ifelse(results$p < 0.05, "Yes", "No")

results$p <- ifelse(
  results$p < 0.001,
  "<0.001",
  sprintf("%.3f", results$p)
)
#sink("t-tests_prot_th.txt", append = TRUE)
results
#sink()

##### NÃO-PARAMETRICO #####

results_w <- data.frame()

for(site in sites){

  test_data <- data[[var]][
    data$Site == site &
    !is.na(data[[var]])
  ]

  wt <- wilcox.test(
    test_data,
    ref,
    alternative = "greater"
  )

  results_w <- rbind(
    results_w,
    data.frame(
      Site = site,
      W = unname(wt$statistic),
      p = wt$p.value
    )
  )
}

results_w$Significant <- ifelse(results_w$p < 0.05, "Yes", "No")

results_w$p <- ifelse(
  results_w$p < 0.001,
  "<0.001",
  sprintf("%.3f", results_w$p)
)

#sink("t-tests_np.txt", append = TRUE)
results_w
#sink()

MXR phenotype activity. Rhodamine B concentrations differed among sampling sites relative to the reference site (Site 1), which exhibited a mean concentration of 38,396.35 ± 4,749.86 RFU mg⁻¹. Significantly higher concentrations were observed at Sites 4, 5 and 6, with mean values of 73,113.69 ± 5,311.65, 67,116.05 ± 7,023.54 and 62,929.41 ± 6,196.21 RFU mg⁻¹, respectively (Wilcoxon rank-sum test; p < 0.01). In contrast, Rhodamine B concentrations at Sites 2 and 3 (16,978.11 ± 2,139.98 and 20,973.12 ± 2,373.55 RFU mg⁻¹, respectively) were actually significantly lower than those observed at the reference site (two-sided Wilcoxon rank-sum test; p < 0.01). The highest mean concentration was recorded at Site 4, corresponding to an approximately 1.9-fold increase relative to the reference site.

Code: Rhodamine B concentration
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$Rodamina_rfu_mg))

##### GRAFICO POR ANÁLISE #####

#N #data <- data %>% mutate(N = 1)
#RPC
#Rodamina_rfu_mg
#TH_musculo_p
#TH_branquias_p
#Proteinas_Totais_mg

# Convert the column to numeric, replace commas with dots for proper conversion
data$Rodamina_rfu_mg <- as.numeric(gsub(",", ".", data$Rodamina_rfu_mg))
VAR <- deparse(substitute(data$Rodamina_rfu_mg))
var <- "Rodamina_rfu_mg" %>% rlang::sym()
str(data)

# Calculate the average for each combination of species, Amostra_tipo, and Site
average <- summarise(
  group_by(data, Site, Especie, Amostra_tipo),
  Average = mean(!!var, na.rm = TRUE),
  n = sum(!is.na(!!var)),
  SE = ifelse(
    sum(!is.na(!!var)) > 1,
    sd(!!var, na.rm = TRUE) / sqrt(sum(!is.na(!!var))),
    0
  ),
  .groups = "drop"
)
average <- average[
  average$Amostra_tipo == "rd" & !is.nan(average$Average),
]
# Plot the average
ggplot(average, aes(x = Especie, y = Average, fill = Amostra_tipo)) +
  geom_bar(stat = "identity", position = position_dodge()) +
  facet_wrap(~ Site) +
  labs(title = paste("Average POR ANÁLISE para", VAR, "per Species and Amostra_tipo for Each Site"),
       x = "Species",
       y = paste("Average", VAR)) +
  theme_minimal() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1))

##### PLOT FINAL #####

rodam <- 
  ggplot(average, aes(x = Especie, y = Average, fill = Amostra_tipo)) +
  geom_bar(stat = "identity", position = position_dodge(), fill = "grey70") +
    geom_errorbar(
    aes(
      ymin = Average - SE,
      ymax = Average + SE
    ),
    width = 0.2
  ) +
  geom_text(
    aes(
      y = 5000,
      label = paste0("n=", n)
    ),
    size = 3
  ) +
  facet_wrap(
  ~ Site,
  labeller = labeller(
    Site = c(
      "B-P05" = "Site 1",
      "B-P06" = "Site 2",
      "B-P07" = "Site 3",
      "B-P08" = "Site 4",
      "B-P09" = "Site 5",
      "B-P10" = "Site 6"
    )
  )
) +
  scale_fill_grey(start = 0.2, end = 0.8) +
  labs(title = "",
       x = "Species",
       y = "Rhodamine B concentration (RFU/mg tissue)",
       fill = "Species") +
  theme_grey() +
  theme(axis.text.x = element_text(angle = 45, hjust = 1),
        legend.position = "none")
rodam
ggsave(rodam, dpi = 300, filename = "fig-rodam.png", bg = "white")
Code: Rhodamine B concentration t-tests: one-sided two-sample t-test
dev.off()
rm(list=ls(all=TRUE))
cat("\014")

data <- read.csv("mxr23_part.csv",
                   sep = ";", dec = ".",
                   row.names = 1,
                   header = TRUE,
                   na.strings = NA)

table(data$Especie, is.na(data$Rodamina_rfu_mg))
var <- "Rodamina_rfu_mg"

# Reference site
ref <- data[[var]][
  data$Site == "B-P05" &
  !is.na(data[[var]])
]

# Sites to compare
sites <- setdiff(unique(data$Site), "B-P05")

# Run tests
results <- data.frame()

for(site in sites){

  test_data <- data[[var]][
  data$Site == site &
  !is.na(data[[var]])
]

  tt <- t.test(
    test_data,
    ref,
    alternative = "greater"
  )

  results <- rbind(
  results,
  data.frame(
    Site = site,
    Mean_Site = mean(test_data),
    SD_Site = sd(test_data),
    SE_Site = sd(test_data) / sqrt(length(test_data)),
    Mean_Ref = mean(ref),
    SD_Ref = sd(ref),
    SE_Ref = sd(ref) / sqrt(length(ref)),
    N_Site = length(test_data),
    N_Ref = length(ref),
    t = unname(tt$statistic),
    df = unname(tt$parameter),
    p = tt$p.value
  )
)
}

results

results$Significant <- ifelse(results$p < 0.05, "Yes", "No")

results$p <- ifelse(
  results$p < 0.001,
  "<0.001",
  sprintf("%.3f", results$p)
)
#sink("t-tests_rodam.txt", append = TRUE)
results
#sink()

##### NÃO-PARAMETRICO #####

results_w <- data.frame()

for(site in sites){

  test_data <- data[[var]][
    data$Site == site &
    !is.na(data[[var]])
  ]

  wt <- wilcox.test(
    test_data,
    ref,
    alternative = "two.sided"
  )

  results_w <- rbind(
    results_w,
    data.frame(
      Site = site,
      W = unname(wt$statistic),
      p = wt$p.value
    )
  )
}

results_w$Significant <- ifelse(results_w$p < 0.05, "Yes", "No")

results_w$p <- ifelse(
  results_w$p < 0.001,
  "<0.001",
  sprintf("%.3f", results_w$p)
)

#sink("t-tests_np.txt", append = TRUE)
results_w
#sink()

4 Discussion

The Atlantic Forest in Brazil has been widely fragmented, and other than larger areas in the southeast of the country, most of it persists as small fragments in the northeast of Brazil (Pinto & Brito, 2003). The importance of these fragments to the aquatic diversity has been reported by a number of authors (Abilhoa et al., 2011; Gouveia et al., 2017; Malabarba, 2006), even though they have often been overlooked (Pereira et al., 2026), along with their ecological services and hydrological health. Our hypothesis that increased pressure and degradation of the surrounding environment of an Atlantic Forest stream can be linked to physiological responses of fish by ecophysiological markers was supported, as shown by a clear increase in all biomarkers from sites outside the protected area, in comparison with the reference site. Furthermore, sites within the protected area did not differ significantly from the reference site, regardless of their type (I or II).

These results demonstrate that, not only the physiological condition of fish can vary consistently with the environmental gradient (Macêdo et al., 2019), but observable changes can be measured in small head water streams, even though they are surrounded by a mosaic of different degrees of management of the riparian zone (Gouveia et al., 2017). Interestingly, Ondei et al. (2020) results show a low spatial response in fish biomarkers, indicating (according to these authors) that fish do not respond to certain spatial scales. Even though they point out that this is a desirable property, since biomonitoring is based on the deviation from a given reference site (Sanchez et al., 2013), in our results, the short scale response of the biomarkers in comparison with the reference site is an advantage, since we were able to detect short-scale changes in fish physiological condition in response to local environmental modifications.

Such modifications were clearly revealed by the environmental variables analyses, which separated between the relatively preserved upstream sites and downstream sites characterized by greater turbidity, increased sandy substrate, reduced overhanging vegetation, and higher concentrations of ammonia and nitrate. It is well registered in the literature the close association between the riparian zone and the stream fauna (Peressin et al., 2020; Terra et al., 2016), yet our study highlights another important association which are changes in the environmental characteristics of the riparian zone and the associated responses in fish physiology to the management and degradation of their surrounding terrestrial habitat in the scale of a small head water stream.

This spatial pattern was also reflected in fish biomarkers, as total plasma protein concentrations and tissue moisture were greater in the downstream sites, which were the most environmentally altered. Interestingly, fish from the reference and other upstream sites (withing the protected area) showed comparatively similar responses (Quiroz-Jara et al., 2021). In a similar study, but using different biomarkers, Cazenave et al. (2009) observed that enven though water quality parametres did differ among study sites, biomarkers responses indicated that fish were living under stressful environmental conditions. On the other hand, Haredi et al. (2020) showed a clear association between low quality of physical and chemical parameters and increased responses of fish biomarkers. In our study, we not only observed measurable physiological responses of fish, but we were also able to discern these responses across the different levels of environmental disturbance. Together, these results suggest that the protection provided by the REBIO is associated with relatively preserved riparian zone, leading to reduced physiological stress, while environmental pressures occurring beyond the protected area are accompanied by detectable responses of the biomarkes evaluated (Ghisi et al., 2022).

Regarding possible temporal patterns, it is worth noting that the present study was conducted during the summer, when inputs from the riparian zone are lower, and stressful conditions are expected to be higher due to the reduced water flow and consequent concentration of chemicals in the water, as also observed by Ondei et al. (2020). During the rainy season, the catchment basin is integrated, and even though there is more runoff from the terrestrial environment, the greater volume of water and connectivity among aquatic habitats may enhance the dilution and transport of dissolved and particulate pollutants (see Sanchez et al. (2007)). Thus, the conditions observed during the summer may represent a period of greater environmental constraint for aquatic organisms, particularly in sites where reduced flow and water retention intensify the effects of local disturbances (Sanchez et al., 2010). This line of thought, as well as, the longitudial effects of unpolluted water from the upper reaches to the more polluted downstream waters, remain to be evaluated, since, mainly during the rainy season, water flow from the upper and proteceted reaches may flush chemical and biological pollutants from the lower reaches.

It was initially expected, in this study, a decrease in muscle total protein and an increase in the tissue water content in fish from disturbed areas (Haredi et al., 2020). The unexpected increase of total plasmatic protein concentration in the animals exposed to pollution, however, is consistent with the Sabae & Mohamed (2015) study, which attributes this increase to five possibilities: activation of metabolic systems as answer to pollutants exposure, degradation of cellular material in the liver, severe pathological conditions in the liver and kidneys, loss of water in the plasma, and induction of proteic synthesis in the liver. Other than the activation of metabolic systems as answer to pollutants, we believe that the liver function may have been impaired in the study fish, as observed by Cazenave et al. (2009). These authors associate oxidative stress induced in the liver to small changes in gills and kidney. It is expected that these different tissues do not respond in the same way to environmental stress, given the higher sensitivity of the liver in comparison with other organs (Cazenave et al., 2006).

The average values of moisture content, parameter related to the osmoregulatory function, were similar for the gill and muscle tissues in the protected area. Regarding the fish from sampling sites subjected to environmental disturbacne, the averages were significantly higher. It is important to note that these increases were only marginal and that these results are not very expressive, like the ones onserved for the other biomakers measured in the present study. Haredi et al. (2020) and Medeiros et al. (2025) observed a similaly low expressiveness in differences in fish muscle moisture content and fiber deterioration between disturbed and non-disturbed sites. Both authors indicate that a range of factors, other than pollution and environmental degradation, will affect fish muscle quality and water content, such as, hardiness of the species, individual variation, maturity state, sex, feeding conditions and the period of fish sampling. Although an important parameter, tissue moisture content is subject to a wide range of factors, making this variable less straightforward to interpret than more precise biomarkers (see Younis et al. (2015)), not being indicated as a parameter to make direct quantitative comparison for animals exposed to environmental degradation, but still useful in making a qualitative evaluation of clearly predefined groups (i.e., fish from impacted versus fish from non-impacted streams).

When we compare the accumulation of rhodamine B concentrations in the gills of the different study species, it is observed different levels of accumulation for some species, such as both Astyanax, indicating different levels of multixenobiotic resistance (MXR) activity for this genera when compared to Hemigrammus living in the same environment. This is an expected result, since the former genus encompasses a more hardy group of species, known to inhabit a hide range of habitats with different levels environmental pressures (Batista et al., 2024; Egea-Serrano et al., 2022; Medeiros et al., 2024). Since the present study does not aim to evaluate species comparatively and all species showed significant differences to Hemigrammus unillineatus found in the reference site, to evaluate the bioindicator potential of the MXR activity in the present study, we compare its activity of all species for each site.

Overall, the present study showed that fish from sites in the area of environmental protection were the ones with the lowest accumulation of rhodamine B. This difference in rhodamine B accumulation between sites within the protected area and sites subject to environmental pollution was opposite to the expected, since the MXR activity is inducible by the presence of xenobiotics in the environment (Smital & Kurelec, 1998).

The expression of this mechanism is mediated mainly by the P‐glycoprotein (P‐gp) (Epel et al., 2008; Macêdo et al., 2019). The P‐gp binds to a variety of substrates and facilitates the efflux of these substrates, thereby preventing the intracellular accumulation of xenobiotics and establishing resistance to these substances (Kurelec, 1995; Minier et al., 1999). Threrefore, the expression of the MXR phenotype could be considered a generalized biomarker of aquatic contamination (Macêdo et al., 2019). In this study, we measured the MXR activity by the use of rhodamine B as a fluorescent substrate of intracellular P‐gp. Because P-gp exports rhodamine B from the cells, high intracellular accumulation indicates reduced P-gp-mediated efflux, and consequently, low MXR activity (Smital & Kurelec, 1998). Meaning that a low accumulation of rhodamine B in the cell is indicative of high MXR activity, whereas a higher accumulation of rhodamine B should indicate a lower MXR activity, which in turn indicates a better environmental condition regarding pollution (Macêdo et al., 2019). This can be explained by the presence of MXR phenotype inhibitors, which were associated with more polluted waters by Kurelec et al. (2000). Among these possible inhibitors Smital et al. (2004) mentions, pesticides, fragrances, microbial degradation products and natural inhibitors from invasive species. These chemosensitizers can cause increase in intracellular accumulation and toxic effects of other xenobiotics normally effluxed by MXR transport proteins (Smital et al., 2004). Since the fish did not go through a decontaminating period, these chemosensitizers could still be bound to the active P-gp site, preventing rhodamine B to attach to it.

Regarding the nutrient concentrations, all variables measured showed very small values. The most informative of these variables were ammonia, nitrate and phosphorus. Even though, present in all sites, ammonia tended to have greater average concentrations and maximum values in sites outsite the conservation unit, which were more prone to receive runoff of agricultural effluentes from the soil. As pointed out by Piedras et al. (2006), elevated concentrations of ammonia (and its nitrification products, such as nitrate) can generate metabolic ammonia retention in fish, causing toxicity. Furthermore, the toxicity of non-ionized ammonia in the aquatic environment is greater for smaller organisms, but also dependent on the interaction with other environmental variables, such as temperature, pH and salinity (Piedras et al., 2006; Randall & Tsui, 2002; Yun et al., 2026). The elevated concentrations of phosphorus in lower stream sites (particularly site 6), is a further indication of contamination of the water by domestic and/or agricultural effluents into the aquatic environment (Johnes et al., 1996; Sabae & Mohamed, 2015; Yan et al., 2016). As a whole the higher concentrations of key nutrients at sites located outside or near the boundaries of the protected area suggest exposure to environmental stressors potentially associated with anthropogenic activities. These results also show an increase of particullar nutrients at lower sites, as pollutants likely accumulate on the sediments and on the riparian zone.

Conclusions

The present study evaluated whether biomarkers of physiological stress in fish from an Atlantic forest head water stream could reflect environmental degradation in the riparian zone. Even though there is no indication of heavy metals in the study basin, pesticides are of current use in the region, indicating the existence of xenobiotics in the water. Our results show that all biomarkers (plasmatic protein concentrations, tissue moisture content and the MXR phenotype activity) responded to the environmental variables indications of environmental degradation and indications of pollution from nutrient concentration. Furthermore, our results indicate that sites within the boundaries of the protected show low to no physiological response from fish that can be associated with environmental stress. Nonetheless, despite a REBIO being one of the most restricted types of protected areas in Brazil, its more immediate vicinity suffers from the historic factors of occupation and absence of an effective buffer zone around the three areas that make up the conservation unit. Our results show that human populations activities cause direct impact to the well being of fish populations, which are a clear reflect of the overall degradation of their surrounding environment. The runoff of contaminants from septic tanks and released agricultural chemicals into the soil from the riparian zone to streams tend to contaminate the surface and groundwater, which led to the observable ecophysiological responses of the ichthyofauna along the gradient of environmental degradation of the study stream. These findings support the use of biochemical and physiological biomarkers as sensitive tools for detecting sublethal effects of environmental degradation, even in landscapes that contain legally protected areas. Finally, in this study we identify physiological biomarkers for fish and associate them with environmental changes across a headwater stream. This will enable decision makers to take action on the influence of anthropogenic threats and the ecological consequence of pollution exposure of the buffer zone.

Acknowledgements

We are grateful to Thainá da Silva Oliveira (Programa de Iniciação Cientifica, PIBIC-UEPB) for assistance with laboratory activities; and to Ellen Gomes da Silva (Programa de Iniciação Cientifica, PIBIC-UEPB) for assistance with fieldwork activities. ACP and ASB are grateful to Programa de Pós-Graduação em Ecologia e Conservação (PPGEC-UEPB) and CAPES for scholarship granted. MMSM and LRCS are grateful to Programa de Iniciação Científica (PIBIC-UEPB) and CNPq for scholarship granted.

Authorship contribution statement

Authorship of this paper is based on CRediT (2026).

Adamastor Coutinho Pinto: Investigation, Formal analysis, Writing – original draft.

Larissa Rafaela Caetano da Silva: Data curation, Formal analysis.

Mayara Mirelly da Silva Monteiro: Data curation, Formal analysis, Software, Writing – review & editing.

Alice da Silva Barros: Data curation, Formal analysis, Software, Writing – review & editing.

Enelise Marcelle Amado: Data curation, Methodology, Resources, Supervision, Validation, Writing – review & editing.

Elvio Sergio Figueredo Medeiros: Conceptualization, Data curation, Formal analysis, Funding acquisition, Methodology, Project administration, Resources, Software, Supervision, Validation, Writing – review & editing.

Ethical approval

Fish collection was authorized by the Instituto Chico Mendes de Conservação da Biodiversidade (license SISBIO-89415-2) and the project was approved by the Ethics Committee on Animal Use of the State University of Paraíba (protocol CEUA-048/2023).

Declaration of generative AI and AI-assisted technologies in the manuscript preparation process

During the preparation of this manuscript, the authors used ChatGPT (OpenAI) to assist with English language editing and improve the clarity and fluency of the writing. No AI was used for data analysis, interpretation, or generation of scientific conclusions. The authors carefully reviewed and revised all AI-generated suggestions and take full responsibility for the accuracy and content of the final manuscript (CNPq, 2026).

References

Abilhoa, V., R., R., Bornatowski, H., & S. Vitule, J. R. (2011). Fishes of the atlantic rain forest streams: Ecological patterns and conservation (O. Grillo, Ed.). InTech. https://doi.org/10.5772/24540
Amoatey, P., & Baawain, M. S. (2019). Effects of pollution on freshwater aquatic organisms. Water Environment Research, 91(10), 1272–1287. https://doi.org/10.1002/wer.1221
APHA. (1998). Standard methods for the examination of water and wastewater. https://www.standardmethods.org/
Arruda, D. B., Cunha, B. P., & Rêgo, K. M. C. (2013). Conflitos entre ReBio guaribas e comunidades locais: (in)justiça ambiental e ecologia política. Revista Direitos Emergentes na Sociedade Global, 2(2).
Ayrapetyan, S. (2012). Cell hydration as a universal marker for detection of environmental pollution. The Environmentalist, 32(2), 210–221. https://doi.org/10.1007/s10669-011-9380-3
Bard, S. (2000). Multixenobiotic resistance as a cellular defense mechanism in aquatic organisms. Aquatic Toxicology, 48(4), 357–389. https://doi.org/10.1016/S0166-445X(00)00088-6
Batista, C. D. P., Ferreira, A. D. A., & Silva, G. J. D. C. (2024). Status of studies about brazilian bioindicator fishes, a review. Acta Limnologica Brasiliensia, 36, e14. https://doi.org/10.1590/s2179-975x8723
Bradford, M. M. (1976). A rapid and sensitive method for the quantitation of microgram quantities of protein utilizing the principle of protein-dye binding.
BRASIL. (2000). Lei no 9.985, de 18 de julho de 2000. Regulamenta o art. 225, § 1o, incisos i, II, III e VII da constituição federal, institui o sistema nacional de unidades de conservação da natureza e dá outras providências. https://www.planalto.gov.br/ccivil_03/leis/l9985.htm
Cazenave, J., Bacchetta, C., Parma, M. J., Scarabotti, P. A., & Wunderlin, D. A. (2009). Multiple biomarkers responses in Prochilodus lineatus allowed assessing changes in the water quality of salado river basin (santa fe, argentina). Environmental Pollution, 157(11), 3025–3033. https://doi.org/10.1016/j.envpol.2009.05.055
Cazenave, J., Bistoni, M. D. L. A., Pesce, S. F., & Wunderlin, D. A. (2006). Differential detoxification and antioxidant response in diverse organs of Corydoras paleatus experimentally exposed to microcystin-RR. Aquatic Toxicology, 76(1), 1–12. https://doi.org/10.1016/j.aquatox.2005.08.011
CNPq. (2026). Portaria 2664/2026 de 6 de março de 2026 - política de integridade na atividade científica do CNPq. http://memoria2.cnpq.br/web/guest/view/-/journal_content/56_INSTANCE_0oED/10157/23142775
Costa, D. A., Azevedo, J. P. S., Santos, M. A., & Assumpção, R. S. F. V. (2020). Water quality assessment based on multivariate statistics and water quality index of a strategic river in the brazilian atlantic forest. Scientific Reports, 10(1), 22038. https://doi.org/10.1038/s41598-020-78563-0
CRediT. (2026). CRediT - contributor role taxonomy. NISO/CRediT Standing Committee. https://credit.niso.org/
David, D. D., Lima, O. G., Nóbrega, A. M. C. D. S., & Amado, E. M. (2018). Capacity of tissue water regulation is impaired in an osmoconformer living in impacted estuaries? Ecotoxicology and Environmental Safety, 166, 375–382. https://doi.org/10.1016/j.ecoenv.2018.09.111
Eerkes-Medrano, D., & Thompson, R. (2018). Occurrence, fate, and effect of microplastics in freshwater systems (E. Y. Zeng, Ed.; pp. 95–132). Elsevier. https://doi.org/10.1016/B978-0-12-813747-5.00004-7
Egea-Serrano, A., Alves, M. C., Solé, M., & Tejedo, M. (2022). Upper thermal tolerances and vulnerability to global warming in a brazilian caatinga fish Astyanax bimaculatus (linnaeus, 1758) population. Austral Ecology, 47(6), 1157–1161. https://doi.org/10.1111/aec.13207
Epel, D. (1998). Use of multidrug transporters as first lines of defense against toxins in aquatic organisms. Comparative Biochemistry and Physiology Part A: Molecular & Integrative Physiology, 120(1), 23–28. https://doi.org/10.1016/S1095-6433(98)10005-3
Epel, D., Luckenbach, T., Stevenson, C. N., MacManus-Spencer, L. A., Hamdoun, A., & Smital, A. T. (2008). Efflux transporters: Newly appreciated roles in protection against pollutants. Environmental Science & Technology, 42(11), 3914–3920. https://doi.org/10.1021/es087187v
Evans, D. H. (1987). The fish gill: Site of action and model for toxic effects of environmental pollutants.
Ghisi, N. D. C., Larentis, C., Oliveira, E. C. D., Neves, M. P., Zavaski, A. G., Roque, A. D. A., Wachtel, C. C., Silva, A. P. D., Lima, E. B. S. D., Costa, G. D. O. N., Amorim, J. P. D. A., Ribeiro, C. A. D. O., & Delariva, R. L. (2022). Environmental assessment of neotropical streams using fish as bioindicators: A multibiomarker and integrated approach. Hydrobiologia, 849(20), 4587–4604. https://doi.org/10.1007/s10750-020-04460-2
Gouveia, R. S. D., Lira, G. L. D. A., Anselmo Ramos, T. P., & Medeiros, E. S. F. (2017). Ichthyofauna of the reserva biológica guaribas and surrounding areas, state of paraíba, brazil. Check List, 13(5), 581–590. https://doi.org/10.15560/13.5.581
Haredi, A. M. M., Mourad, M., Tanekhy, M., Wassif, E., & Abdel-Tawab, H. S. (2020). Lake edku pollutants induced biochemical and histopathological alterations in muscle tissues of nile tilapia (oreochromis niloticus). Toxicology and Environmental Health Sciences, 12(3), 247–255. https://doi.org/10.1007/s13530-020-00042-w
Harmon, S. M. (2018). The effects of microplastic pollution on aquatic organisms (E. Y. Zeng, Ed.; pp. 249–270). Elsevier. https://doi.org/10.1016/B978-0-12-813747-5.00008-4
Heinrichs, R., Otto, R., Magalhães, A., & Meirelles, G. C. (2017). Importance of sugarcane in brazilian and world bioeconomy (S. Dabbert, I. Lewandowski, J. Weiss, & A. Pyka, Eds.; pp. 205–217). Springer International Publishing. https://doi.org/10.1007/978-3-319-58374-7_11
Johnes, P., Moss, B., & Phillips, G. (1996). The determination of total nitrogen and total phosphorus concentrations in freshwaters from land use, stock headage and population data: Testing of a model for use in conservation and water quality management. Freshwater Biology, 36(2), 451–473. https://doi.org/10.1046/j.1365-2427.1996.00099.x
Kurelec, B. (1992). The multixenobiotic resistance mechanism in aquatic organisms. Critical Reviews in Toxicology, 22(1), 23–43. https://doi.org/10.3109/10408449209145320
Kurelec, B. (1995). Inhibition of multixenobiotic resistance mechanism in aquatic organisms: Ecotoxic consequences. Science of The Total Environment, 171(1-3), 197–204. https://doi.org/10.1016/0048-9697(95)04689-4
Kurelec, B., Smital, T., Pivèeviæ, B., Eufemia, N., & Epel, D. (2000). Multixenobiotic resistance, p-glycoprotein, and chemosensitizers. Ecotoxicology, 9(5), 307–327. https://doi.org/10.1023/A:1026560922731
Lê, S., Josse, J., & Husson, F. (2008). FactoMineR: An r package for multivariate analysis. Journal of Statistical Software, 25(1), 1–18.
Legendre, P., & Legendre, L. (1998). Numerical ecology (2nd English ed).
Macêdo, A. K. S., Da Silva, J. R. P., Dos Santos, H. B., Thomé, R. G., Vendel, A. L., & Amado, E. M. (2019). Estuarine fish assemblages present a species-specific difference in the multixenobiotics resistance activity. Journal of Experimental Zoology Part A: Ecological and Integrative Physiology, 331(10), 530–539. https://doi.org/10.1002/jez.2320
Mackereth, F. J. H., Heron, J., & Talling, J. F. (1978). Water analysis: Some revised methods for limnologists.(scientific publication n. 36) freshwater biological association. Cumbria and Dorset, England.
Maitland, P. S. (1990). Field studies: Sampling in freshwaters (pp. 123–148). Springer US. https://doi.org/10.1007/978-1-4613-0501-9_5
Malabarba, L. R. (2006). Peixes de riachos da mata atlântica nas unidades de conservação do vale do rio ribeira de iguape no estado de são paulo. Neotropical Ichthyology, 4(3), 377–377. https://doi.org/10.1590/S1679-62252006000300012
Marques, M. C. M., & Grelle, C. E. V. (Eds.). (2021). The atlantic forest: History, biodiversity, threats and opportunities of the mega-diverse forest. Springer International Publishing. https://doi.org/10.1007/978-3-030-55322-7
Marques, M. C. M., Trindade, W., Bohn, A., & Grelle, C. E. V. (2021). The atlantic forest: An introduction to the megadiverse forest of south america (M. C. M. Marques & C. E. V. Grelle, Eds.; pp. 3–23). Springer International Publishing. https://doi.org/10.1007/978-3-030-55322-7_1
McCune, B., & Grace, J. B. (2002). Analysis of ecological communities. MjM Software Design.
Medeiros, E. S. F., Monteiro, M. M. D. S., Lustosa-Costa, S. Y., Falkenberg, J. M., & Anselmo Ramos, T. P. (2025). Development of lip protrusion in trahira, hoplias malabaricus (bloch, 1794), associated with aquatic surface respiration. Neotropical Biology and Conservation, 20(2), 103–113. https://doi.org/10.3897/neotropical.20.e146972
Medeiros, E. S. F., Silva, M. J. D., Ramos, T. P. A., & Ramos, R. T. C. (2024). Environmental variables as predictors of fish community composition in semiarid aquatic systems. Acta Limnologica Brasiliensia, 36, e4. https://doi.org/10.1590/s2179-975x3023
Medeiros, E. S. F., Silva, M., Figueiredo, B., Ramos, T., & Ramos, R. (2010). Effects of fishing technique on assessing species composition in aquatic systems in semi-arid brazil. Brazilian Journal of Biology, 70, 255–262. https://doi.org/10.1590/S1519-69842010000200004
Medeiros, E. S. F., Silva, M. J., & Ramos, R. T. C. (2008). Application of catchment- and local-scale variables for aquatic habitat characterization and assessment in the brazilian semi-arid region. Neotropical Biology and Conservation, 3(1), 13–20.
Minier, C., Eufemia, N., & Epel, D. (1999). The multi-xenobiotic resistance phenotype as a tool to biomonitor the environment. Biomarkers, 4(6), 442–454. https://doi.org/10.1080/135475099230606
Miranda, J. (2012). Ameaças aos peixes de riachos da mata atlântica. Ameaças Aos Peixes de Riachos Da Mata Atlântica.
MMA. (2003). Plano de manejo reserva biológica guaribas. Fase 2. https://bit.ly/4gtUh3i
MMA. (2026). Painel de unidades de conservação brasileras. https://encurtador.com.br/kLle
Mugodo, J., Kennard, M. J., Liston, P., Nichols, S., Linke, S., Norris, R. H., & Lintermans, M. (2006). Local stream habitat variables predicted from catchment scale characteristics are useful for predicting fish distribution. Hydrobiologia, 572(1), 59–70. https://doi.org/10.1007/s10750-006-0252-7
Myers, N., Mittermeier, R. A., Mittermeier, C. G., Da Fonseca, G. A. B., & Kent, J. (2000). Biodiversity hotspots for conservation priorities. Nature, 403(6772), 853–858. https://doi.org/10.1038/35002501
Ondei, L. D. S., Teresa, F. B., Garcia, D. P., Felício, A. A., Da Silva, D. G. H., & De Almeida, E. A. (2020). Fish biomarker responses to perturbation by drought in streams. Neotropical Ichthyology, 18(2), e190130. https://doi.org/10.1590/1982-0224-2019-0130
Padial, A. A., Costa, A. P. L., Bonecker, C. C., Nogueira, D. G., Roque, F. O., Message, H. J., Diniz, L. P., Cavalcante, L. L., Vallejos, M. A. V., Ota, R. R., & Ludwig, T. A. V. (2021). Freshwater studies in the atlantic forest: General overview and prospects (M. C. M. Marques & C. E. V. Grelle, Eds.; pp. 205–230). Springer International Publishing. https://doi.org/10.1007/978-3-030-55322-7_10
Peel, M. C., Finlayson, B. L., & McMahon, T. A. (2007). Updated world map of the köppen-geiger climate classification. Hydrol. Earth Syst. Sci.
Pereira, K. M. G., Cordeiro, N. G., Suzuki, L. P. Z. L., & Berg, E. van den. (2026). How important are small fragments for restoration in the central portion of the atlantic forest? Biological Conservation, 321, Article 111958. https://doi.org/https://doi.org/10.1016/j.biocon.2026.111958
Peressin, A., Casarim, R., Prado, I. G., & Cetra, M. (2020). Physical habitat as predictor of fish trophic structure in brazilian atlantic rainforest streams. Neotropical Ichthyology, 18(2), e190076. https://doi.org/10.1590/1982-0224-2019-0076
Piedras, S., Oliveira, J., Moraes, P., & Bager, A. (2006). Toxicidade aguda da amônia não ionizada e do nitrito em alevinos de cichlasoma facetum (jenyns, 1842). Ciência e Agrotecnologia, 30, 1008–1012. https://doi.org/10.1590/S1413-70542006000500027
Pinto, L. P., & Brito, M. C. W. (2003). Dynamics of biodiversity loss in the brazilin atlantic forest: An introduction (C. Galindo Leal & I. G. Câmara, Eds.; pp. 27–30). Island Press.
Pusey, B., Kennard, M. J., & Arthington, A. (2004). Study area, data collection, analysis and presentation (B. Pusey, M. J. Kennard, & A. Arthington, Eds.; pp. 26–48). CSIRO Publishing.
Quiroz-Jara, M., Casini, S., Fossi, M. C., Orrego, R., Gavilán, J. F., & Barra, R. (2021). Integrated physiological biomarkers responses in wild fish exposed to the anthropogenic gradient in the biobío river, south-central chile. Environmental Management, 67(6), 1145–1157. https://doi.org/10.1007/s00267-021-01465-y
R Core Team. (2017). R: A language and environment for statistical computing. https://www.R-project.org/
R Studio Team. (2022). RStudio: Integrated development environment for r. RStudio, PBC. https://posit.co/products/open-source/rstudio/
Randall, D. J., & Tsui, T. K. N. (2002). Ammonia toxicity in fish. Marine Pollution Bulletin, 45(1-12), 17–23. https://doi.org/10.1016/S0025-326X(02)00227-8
Ross, A., & Willson, V. L. (2017). Basic and advanced statistical tests. SensePublishers. https://doi.org/10.1007/978-94-6351-086-8
Sabae, S. Z., & Mohamed, F. A. S. (2015). Effect of environmental pollution on the health of tilapia spp. From lake qarun.
Sanchez, W., Aït-Aïssa, S., Palluel, O., Ditche, J.-M., & Porcher, J.-M. (2007). Preliminary investigation of multi-biomarker responses in three-spined stickleback (gasterosteus aculeatus l.) sampled in contaminated streams. Ecotoxicology, 16(2), 279–287. https://doi.org/10.1007/s10646-006-0131-z
Sanchez, W., Burgeot, T., & Porcher, J.-M. (2013). A novel integrated biomarker response calculation based on reference deviation concept. Environmental Science and Pollution Research, 20(5), 2721–2725. https://doi.org/10.1007/s11356-012-1359-1
Sanchez, W., Piccini, B., Maillot-Maréchal, E., & Porcher, J.-M. (2010). Comparison of two reference systems for biomarker data analysis in a freshwater biomonitoring context. Environment International, 36(4), 377–382. https://doi.org/10.1016/j.envint.2010.02.006
Santos, M. B., Monteiro Neto, I. E., Souza Melo, S. R. C., & Amado, E. M. (2017). Hemolymph and gill carbonic anhydrase are more sensitive to aquatic contamination than mantle carbonic anhydrase in the mangrove oyster crassostrea rhizophorae. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology, 201, 19–25. https://doi.org/10.1016/j.cbpc.2017.08.008
Smital, T., & Kurelec, B. (1998). The chemosensitizers of multixenobiotic resistance mechanism in aquatic invertebrates: A new class of pollutants. Mutation Research - Fundamental and Molecular Mechanisms of Mutagenesis, 399(1), 43–53. https://doi.org/10.1016/S0027-5107(97)00265-0
Smital, T., Luckenbach, T., Sauerborn Klobucar, R., Hamdoun, A., Vega, R., & Epel, D. (2004). Emerging contaminants - pesticides, PPCPs, microbial degradation products and natural substances as inhibitors of multixenobiotic defense in aquatic organisms. Mutation Research, 552, 101–117. https://doi.org/10.1016/j.mrfmmm.2004.06.006
Soler, J. M. P. (2004). Planejamento de experimentos e pesquisa em limnologia (pp. 16–24). BICUDO, C. E. M. BICUDO, D. C.
Terra, B. D. F., Hughes, R. M., & Araújo, F. G. (2016). Fish assemblages in atlantic forest streams: The relative influence of local and catchment environments on taxonomic and functional species. Ecology of Freshwater Fish, 25(4), 527–544. https://doi.org/10.1111/eff.12231
Yan, Z., Han, W., Peñuelas, J., Sardans, J., Elser, J. J., Du, E., Reich, P. B., & Fang, J. (2016). Phosphorus accumulates faster than nitrogen globally in freshwater ecosystems under anthropogenic impacts. Ecology Letters, 19(10), 1237–1246. https://doi.org/10.1111/ele.12658
Younis, E. M., Al-Asgah, N. A., Abdel-Warith, A.-W. A., & Al-Mutairi, A. A. (2015). Seasonal variations in the body composition and bioaccumulation of heavy metals in nile tilapia collected from drainage canals in al-ahsa, saudi arabia. Saudi Journal of Biological Sciences, 22(4), 443–447. https://doi.org/10.1016/j.sjbs.2014.11.020
Yun, S. C., Jeong, H., Lee, J.-S., Kim, J.-H., Kim, I.-C., Maszczyk, P., Yang, Z., Hagiwara, A., & Lee, J.-S. (2026). A review of ammonia toxicity on aquatic organisms: Species-specific responses, microbial shifts, and environmental interactions. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology, 300, 110388. https://doi.org/10.1016/j.cbpc.2025.110388
Žaja, R., Klobučar, G. I. V., Sauerborn Klobučar, R., Hackenberger, B. K., & Smital, T. (2006). Haemolymph as compartment for efficient and non-destructive determination of p-glycoprotein (pgp) mediated MXR activity in bivalves. Comparative Biochemistry and Physiology Part C: Toxicology & Pharmacology, 143(1), 103–112. https://doi.org/10.1016/j.cbpc.2005.12.009
Zar, J. H. (2010). Biostatistical analysis (5th ed.). Prentice Hall.

Figures and Tables

Figure 1: Location of sampling sites along the Barro Branco Stream (Mamanguape, PB) (Site 1 06°43’06”S and 35°10’54”W; Site 2 06°42’38”S and 35°10’38”W; Site 3 06°43’06”S and 35°10’54”W; Site 4 06°41’46”S and 35°10’36”W; Site 5 06°40’56”S and 35°10’27”; Site 6 06°40’18”S and 35°10’34”W). *Site 1 is the reference site for this study.

Figure 2: Images of the main four types of sites found along the Barro Branco stream (Mamanguape, PB). Type I, well preserved and within the conservation unit. Type II site, preserved and outside the conservation unit. Type III site, disturbed but withing the conservation unit. Type IV site, disturbed and outside the conservation unit.

Figure 3: Principal Component Analysis of the environmental variables for the Barro Branco stream (Mamanguape, PB).

Figure 4: Counts of number of individuals per species and biomarker use (MXR, Multixenobiotic resistance phenotype activity; Moisture, Plasmatic protein dosage and moisture content) for each study site.

Figure 5: Total plasmatic proteins (mg) variance by sampling site and fish species. An * indicates a significant difference from the reference site (1), Wilcoxon rank-sum testa at p < 0.05.

Figure 6: Gill and muscle moisture content (%) by sampling site and fish species. An * indicates a significant difference from the reference site (1), Wilcoxon rank-sum testa at p < 0.05.

Figure 7: Average fluorescence of rhodamin B expressed by relative fluorescence units/milligrams of gill tissue of fish by sampling site. An * indicates a significant difference from the reference site (1), Wilcoxon rank-sum testa at p < 0.05.

Table 1: Environmental variables averaged for sampling sites (minimum-maximum) for the Barro Branco stream (Mamanguape, PB). Site 1 is the reference site for this study.

Table 2: Nutrient concentration (μg L⁻¹) averaged for sampling sites (minimum-maximum) for the Barro Branco stream (Mamanguape, PB). Site 1 is the reference site for this study.

Figures

Fig.: Location of sampling sites
Figure 1: Location of sampling sites.
Fig.: Images of the main four types of sites
Figure 2: Images of the main four types of sites.
Fig.: Principal Component Analysis of the environmental variables
Figure 3: Principal Component Analysis of the environmental variables.
Fig.: Counts of number of individuals per species and biomarker
Figure 4: Counts of number of individuals per species and biomarker.
Fig.: Total protein variance
Figure 5: Total plasmatic proteins variance.
Fig.: Moisture content of tissues
Figure 6: Gill and muscle moisture content.
Fig.: Average fluorescence for rhodamin B
Figure 7: Average fluorescence of rhodamin B.

Tables

Tab.: Environmental data table
Table 1: Environment data GT table.
Site Site 1* Site 2 Site 3 Site 4 Site 5 Site 6
Type I I II III IV IV
Water quality





Temperature (°C) 25.8 (24.9-27.3) 26.5 (25.3-27.5) 25.2 (24.9-25.5) 27.8 (27-28.5) 26.9 (26.6-27.3) 26.3 (26.2-26.5)
Dissolved oxygen (mg/L) 8.1 (7.8-8.5) 7.5 (7.2-7.8) 7.9 (7.7-8) 7.5 (6.5-8.4) 7.9 (7.8-7.9) 5.4 (4.9-6.1)
Conductivity (uS/cm) 146.4 (134-164.4) 148.3 (134.1-156.8) 129 (98.3-160.2) 160.6 (158.8-161.8) 163.3 (161.2-165.5) 154.3 (149.7-157.9)
Turbidity (NTU) 7.4 (4.6-10.7) 6.1 (2.1-14.1) 4.4 (1.5-9.2) 6.4 (6.1-6.8) 17.5 (17.2-17.7) 29.6 (24.2-34.8)
pH 7 (6.3-7.5) 6.8 (5.9-7.7) 5.5 (4.6-6) 5.9 (5.9-6) 6.6 (6.6-6.6) 6.9 (6.6-7.3)
Morphology





Elevation (m a.s.l.) 107.1 (107.1-107.1) 112.8 (112.8-112.8) 89 (89-89) 89.1 (89.1-89.1) 71 (71-71) 41.5 (41.5-41.5)
Width (m) 3.6 (3.6-3.6) 3.7 (2.3-4.4) 1.7 (1.1-2.1) 3.7 (2-4.6) 1.2 (0.8-1.6) 7.2 (7.2-7.2)
Water velocity (m/s) 0.025 (0.025-0.025) 0.035 (0.035-0.035) 0.073 (0.022-0.108) 0.038 (0-0.115) 0.293 (0.25-0.348) 0.022 (0.022-0.022)
Slope 1.2 (1-1.5) 2.5 (2.5-2.5) 2.2 (1.5-2.5) 2 (2-2) 2.7 (2-3) 1.5 (1-2)
Average depth (cm) 13.7 (13.7-13.7) 23 (20.3-28.3) 19.7 (16-22.7) 44.2 (30-51.3) 16.8 (15-19.7) 40.7 (40.7-40.7)
Substrate composition (%)





Mud 100 (100-100) 100 (100-100) 90 (90-90) 80 (80-80) 70 (60-80) 30 (20-50)
Sand 0 0 10 (10-10) 20 (20-20) 30 (20-40) 70 (50-80)
Habitat structures (%)





Aquatic macrophytes 5 (5-5) 20 (0-30) 0 35.8 (22.5-60) 0 1.7 (0-5)
Littoral grass 0 5 (5-5) 2 (0-5) 39.2 (20-50) 0 13.7 (1-30)
Submerged vegetation 5 (5-5) 8.3 (5-15) 3.3 (0-10) 3.3 (0-5) 16.7 (10-20) 5 (0-10)
Overhanging vegetation 90 (90-90) 81.7 (80-85) 85 (80-90) 30 (30-30) 70 (60-80) 8.3 (0-20)
Leaf litter 10 (10-10) 3.7 (1-5) 20 (5-30) 26.7 (20-30) 13.7 (1-30) 5.3 (1-10)
Algae 5 (5-5) 0.3 (0-1) 11 (2-21) 16.7 (15-20) 0 (0-0) 0 (0-0)
Debris 30 (30-30) 5.3 (5-6) 10.3 (5-15) 31.7 (30-35) 17 (1-40) 15 (10-20)
Root masses 0 8.3 (5-15) 10 (5-15) 0 76.7 (50-90) 3.3 (0-10)
Tab.: Nutrient data table
Table 2: Nutrient concentrations at the sampling sites.
Variable Site 1 Site 2 Site 3 Site 4 Site 5 Site 6
Ammonia (µg L⁻¹) 54.773 (39.32-64.32) 58.867 (47.5-77.05) 72.047 (33.41-115.68) 81.587 (68.86-92.95) 80.833 (65.23-99.32) 243.71 (203.86-277.5)
Nitrite (µg L⁻¹) 0 (0-0) 0 (0-0) 0 (0-0) 0 (0-0) 0 (0-0) 0 (0-0)
Nitrate (µg L⁻¹) 6 (0-9.7) 8.8 (6.9-12.6) 5.5 (0-12.6) 18.8 (11.1-31.1) 18.3 (16.9-19.7) 10.2 (5.4-16.9)
Orthophosphate (µg L⁻¹) 0 (0-0) 0 (0-0) 1 (0-3) 0 (0-0) 0 (0-0) 0 (0-0)
Phosphorus (µg L⁻¹) 0 (0-0) 4.6 (0.3-10.3) 50.1 (0-150.3) 9 (0-27) 19.2 (15.3-22) 38.7 (13.7-55.3)

Apendices

Non-used figures and tables

WarningNon-used codes
non-used code
table(1:10)
Code: Environment data GT table
library(readr)
library(dplyr)
library(gt)

m_wide_hab <- read.table("m_wide_hab.txt",
                  sep = ";", dec = ".",
                  row.names = 1,
                  header = TRUE,
                  na.strings = NA)

#m_mapa <- read_tsv("column_labels.txt", show_col_types = FALSE)
dados <- m_wide_hab
dados
nomes <- names(dados)[-1]
nomes

df_nomes <- data.frame(
  original = nomes,
  final    = nomes,
  stringsAsFactors = FALSE)

df_nomes
df_nomes$final <- gsub("_", "<br>", df_nomes$final)
dados <- mutate(dados,across(-Variable, ~ gsub(" \\(", "<br>(", .)))
#fix(df_nomes)
#write.table(df_nomes, "df_nomes.txt")
df_nomes <- read.table("df_nomes.txt")

labels_finais <- setNames(
  lapply(df_nomes$final, md),
  nomes)

tabela_gt <- dados %>%
  gt(rowname_col = "Variable") %>%
#  tab_header(
#    title = "Características Ambientais",
#    subtitle = "Valores apresentados como média (mín–máx)"
#  ) %>%
  fmt_markdown(
    columns = -Variable) %>%
  cols_align(
    align = "right",
    columns = -Variable) %>%
  cols_label(.list = labels_finais)

tabela_gt
gtsave(tabela_gt, "gt-hab_gttable.html")
gtsave(tabela_gt, "fig-hab_gttable.png")
saveRDS(tabela_gt, "gt-hab_gttable.rds")