diff --git a/environment-dev.yml b/environment-dev.yml
index 9511ef47..2a14ba65 100644
--- a/environment-dev.yml
+++ b/environment-dev.yml
@@ -3,10 +3,9 @@ channels:
   - conda-forge
 dependencies:
   - python >=3.9,<3.12
-  # Don't forget to sync the changes here with environment.yml!
+  # Don't forget to sync changes between environment.yml, environment-dev.yml, and setup.py!
   # Main packages
   - numpy
-  - scipy
   - statsmodels
   - xarray
   - xclim >=0.45.0
diff --git a/environment.yml b/environment.yml
index 98b2726c..c088fd45 100644
--- a/environment.yml
+++ b/environment.yml
@@ -2,9 +2,9 @@ name: xhydro
 channels:
   - conda-forge
 dependencies:
+  # Don't forget to sync changes between environment.yml, environment-dev.yml, and setup.py!
   - python >=3.9,<3.12
   - numpy
-  - scipy
   - statsmodels
   - xarray
   - xclim >=0.45.0
diff --git a/setup.py b/setup.py
index a46e14d0..94dc808f 100644
--- a/setup.py
+++ b/setup.py
@@ -7,7 +7,8 @@
 with open("README.rst") as readme_file:
     readme = readme_file.read()
 
-requirements = ["numpy", "scipy", "statsmodels", "xarray", "xclim>=0.45.0", "xscen"]
+# Don't forget to sync changes between environment.yml, environment-dev.yml, and setup.py!
+requirements = ["numpy", "statsmodels", "xarray", "xclim>=0.45.0", "xscen"]
 
 dev_requirements = ["pytest", "pytest-cov"]
 
diff --git a/xhydro/frequency_analysis/local.py b/xhydro/frequency_analysis/local.py
index a0017aaf..db32efea 100644
--- a/xhydro/frequency_analysis/local.py
+++ b/xhydro/frequency_analysis/local.py
@@ -1,10 +1,9 @@
 """Local frequency analysis functions and utilities."""
 
 import datetime
-from typing import Union
+from typing import Optional, Union
 
 import numpy as np
-import scipy.stats
 import statsmodels
 import xarray as xr
 import xclim.indices.stats
@@ -18,7 +17,10 @@
 
 
 def fit(
-    ds, distributions: list = None, min_years=None, method: str = "ML"
+    ds,
+    distributions: Optional[list[str]] = None,
+    min_years: Optional[int] = None,
+    method: str = "ML",
 ) -> xr.Dataset:
     """Fit multiple distributions to data.
 
@@ -26,10 +28,10 @@ def fit(
     ----------
     ds : xr.Dataset
         Dataset containing the data to fit. All variables will be fitted.
-    distributions : list of str
+    distributions : list of str, optional
         List of distribution names as defined in `scipy.stats`. See https://docs.scipy.org/doc/scipy/reference/stats.html#continuous-distributions.
         Defaults to ["expon", "gamma", "genextreme", "genpareto", "gumbel_r", "pearson3", "weibull_min"].
-    min_years : int
+    min_years : int, optional
         Minimum number of years required for a distribution to be fitted.
     method : str
         Fitting method. Defaults to "ML" (maximum likelihood).
@@ -86,7 +88,9 @@ def fit(
     return out
 
 
-def parametric_quantiles(p, t: Union[float, list], mode: str = "max") -> xr.Dataset:
+def parametric_quantiles(
+    p: xr.Dataset, t: Union[float, list[float]], mode: str = "max"
+) -> xr.Dataset:
     """Compute quantiles from fitted distributions.
 
     Parameters
@@ -159,7 +163,7 @@ def parametric_quantiles(p, t: Union[float, list], mode: str = "max") -> xr.Data
     return out
 
 
-def criteria(ds, p) -> xr.Dataset:
+def criteria(ds: xr.Dataset, p: xr.Dataset) -> xr.Dataset:
     """Compute information criteria (AIC, BIC, AICC) from fitted distributions, using the log-likelihood.
 
     Parameters
diff --git a/xhydro/indicators.py b/xhydro/indicators.py
index 52f51e67..9123e8a1 100644
--- a/xhydro/indicators.py
+++ b/xhydro/indicators.py
@@ -1,5 +1,6 @@
 """Module to compute indicators using xclim's build_indicator_module_from_yaml."""
 import warnings
+from typing import Optional
 
 import xarray as xr
 import xclim as xc
@@ -17,7 +18,7 @@
 
 
 def compute_volume(
-    da: xr.DataArray, *, out_units: str = "m3", attrs: dict = None
+    da: xr.DataArray, *, out_units: str = "m3", attrs: Optional[dict] = None
 ) -> xr.DataArray:
     """Compute the volume of water from a streamflow variable, keeping the same frequency.
 
@@ -27,7 +28,7 @@ def compute_volume(
         Streamflow variable.
     out_units : str
         Output units. Defaults to "m3".
-    attrs : dict
+    attrs : dict, optional
         Attributes to add to the output variable.
         Default attributes for "long_name", "units", "cell_methods" and "description" will be added if not provided.
 
@@ -58,9 +59,9 @@ def get_yearly_op(
     *,
     input_var: str = "streamflow",
     window: int = 1,
-    timeargs: dict = None,
+    timeargs: Optional[dict] = None,
     missing: str = "skip",
-    missing_options: dict = None,
+    missing_options: Optional[dict] = None,
     interpolate_na: bool = False,
 ) -> xr.Dataset:
     """
@@ -77,7 +78,7 @@ def get_yearly_op(
     window: int
         Size of the rolling window. A "mean" operation is performed on the rolling window before the call to xclim.
         This parameter cannot be used with the "sum" operation.
-    timeargs: dict
+    timeargs: dict, optional
         Dictionary of time arguments for the operation. Keys are the name of the period that will be added to the results (e.g. "winter", "summer", "annual").
         Values are up to two dictionaries, with both being optional.
         The first is {'freq': str}, where str is a frequency supported by xarray (e.g. "YS", "AS-JAN", "AS-DEC"). It needs to be a yearly frequency. Defaults to "AS-JAN".
@@ -87,7 +88,7 @@ def get_yearly_op(
     missing: str
         How to handle missing values. One of "skip", "any", "at_least_n", "pct", "wmo".
         See :py:func:`xclim.core.missing` for more information.
-    missing_options: dict
+    missing_options: dict, optional
         Dictionary of options for the missing values' method. See :py:func:`xclim.core.missing` for more information.
     interpolate_na: bool
         Whether to interpolate missing values before computing the operation. Only used with the "sum" operation. Defaults to False.