Can you help me understand how to fix the ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments

[ben29med]

I had a problem with L-moments when i excute my climate indices command: rocess_climate_indices --index spi --periodicity monthly --netcdf_precip sliced_pre_1901_18.nc --var_name_precip pre --output_file_base sliced_pr.nc --scales 12 --calibration_start_year 1901 --calibration_end_year 2018 What's the problem ! @monocongo @monocongo @WeatherGod @ScottWales @dawiedotcom

[monocongo]

Please post a link to your datasets so I can attempt to reproduce and understand the error. For example sliced_pre_1901_18.nc There have been many issues regarding L-moments in the past, and all seem to have been resolved now if using the latest scipy etc. Please install the climate_indices package from the master branch (which is ahead of what's now on PyPI) to see if that will fix your issue. My apologies, I haven't had time to update PyPI with the latest/greatest code yet.

[ben29med]

I have now an output for spi but with masked array, when I visualize the nc file there are a warning in my data masked array (warning valid max and min). I'll send you my output spi to check if my output is correct or not. Are you recieve my output drive link ?

On Mon, May 22, 2023, 15:58 James Adams @.***> wrote:

Please post a link to your datasets so I can attempt to reproduce and understand the error. For example sliced_pre_1901_18.nc There have been many issues regarding L-moments in the past, and all seem to have been resolved now if using the latest scipy etc. Please install the climate_indices package from the master branch (which is ahead of what's now on PyPI) to see if that will fix your issue. My apologies, I haven't had time to update PyPI with the latest/greatest code yet.

— Reply to this email directly, view it on GitHub https://github.com/monocongo/climate_indices/issues/512#issuecomment-1557373087, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJ5OUUYMEMIQXI6JBVTJS7LXHN5INANCNFSM6AAAAAAXVXSRYU . You are receiving this because you authored the thread.Message ID: @.***>

[monocongo]

Please post a link to your data here.

Are you using the latest code from the master branch of this repository, or are you using the version from PyPI which is a little behind that (sorry for out of sync, tech debt)? Please use the latest code from this repository which should have the L-moments issue worked out. There have been many problems with the L-moments calculations over the years but seems like they've been rectified now.

[caiomattos]

I'm having a similar issue, which I've managed to trace to a violation of theorem 1 in Hosking (1990), which bounds tau3 (or t3 in the code here) to -1 < t3 < 1. In my case, I get t3 = 1. Below is a summary of how I found it.

I have precipitation and evapotranspiration data with shape (lat, lon, time) = (8520, 7320, 180), in monthly time steps from Jan 2004 to Dec 2018. I'm trying to calculate SPEI.

I added some print messages to the main loop in main._apply_along_axis_double() to get the lat, lon indexes of where the L-moments error were happening. Then I proceeded with a step-by-step calculation using the modules. Below is an example of what is causing the problem.

For the month of August, we have:

precipitation = [0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0., 0.]
evapotranspiration = [1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1., 1.]
time_step_values = [2997., 2997., 2997., 2997., 2997., 2994., 2997., 2997., 2997., 2997., 2997., 2997., 2997., 2997., 2997.]
lmoments = [ 2.9968e+03  2.0000e-01 -1.0000e+00]

When the calculation calls lmoments._estimate_pearson3_parameters(), we get

t3 = abs(lmoments[2])  # L-skewness?
t3 >= 1
True

And hence, the error. Using scipy.stats.pearson3.fit() on the data with both maximum likelihood or method of moments (MM) we get valid output:

scipy.stats.pearson3.fit(time_step_values, method='MM')
(-3.474284756042035, 2996.800000000001, 0.7483314723274126)

Just for comparison, for the month of February, which doesn't trigger an error, we get estimates that are pretty close between scipy with "MM" and the estimate_pearson3_parameters routine:

scipy: {'loc': 3008.3333333333703, 'skew': 0.5728926162569001, 'scale': 4.422166348050119}
estimate_pearson3_parameters: {'loc': 3008.3333333333335, 'skew': 0.6312923338011005, 'scale': 4.68319994913074}

Looking at other points this seems to happen when all P = 0 for a certain month (where seasonality is high). Maximum likelihood has been used before to calculate SPI/SPEI, with L-moments providing initial values (as in Stagge et al. 2015, https://doi.org/10.1002/joc.4267), so this might be an option in the future for similar cases.

[monocongo]

@caiomattos' analysis looks really solid here. Thanks for the input! If there's a way to fix our code to accommodate this then a PR will be welcome.

[Darren-Ray]

Hi,

I installed climate-indices into a fresh Python3.10 virtual env with the latest scipy version 1.10.1 and climate indices 1.10.13 and am getting the same error: Traceback (most recent call last): File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 125, in worker result = (True, func(*args, *kwds)) File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 48, in mapstar return list(map(args)) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1428, in _apply_along_axis_double computed_array[i, j] = func1d(x[j], y[j], parameters=params["args"]) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1151, in _spei return indices.spei(precips_mm=precips, File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/indices.py", line 333, in spei compute.transform_fitted_pearson( File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/compute.py", line 553, in transform_fitted_pearson pearson_parameters( File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/compute.py", line 349, in pearson_parameters params = lmoments.fit(time_step_values) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/lmoments.py", line 30, in fit raise ValueError(message) ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments """

The above exception was the direct cause of the following exception:

Traceback (most recent call last): File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1782, in main _compute_write_index(kwrgs) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1044, in _compute_write_index _parallel_process(keyword_arguments["index"], File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1324, in _parallel_process pool.map(_apply_along_axis_double, chunk_params) File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 364, in map return self._map_async(func, iterable, mapstar, chunksize).get() File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 771, in get raise self._value ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments multiprocessing.pool.RemoteTraceback: """ Traceback (most recent call last): File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 125, in worker result = (True, func(*args, *kwds)) File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 48, in mapstar return list(map(args)) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1428, in _apply_along_axis_double computed_array[i, j] = func1d(x[j], y[j], parameters=params["args"]) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1151, in _spei return indices.spei(precips_mm=precips, File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/indices.py", line 333, in spei compute.transform_fitted_pearson( File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/compute.py", line 553, in transform_fitted_pearson pearson_parameters( File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/compute.py", line 349, in pearson_parameters params = lmoments.fit(time_step_values) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/lmoments.py", line 30, in fit raise ValueError(message) ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments """

The above exception was the direct cause of the following exception:

Traceback (most recent call last): File "/home/dray/anaconda3/envs/CI/bin/process_climate_indices", line 8, in sys.exit(main()) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1782, in main _compute_write_index(kwrgs) File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1044, in _compute_write_index _parallel_process(keyword_arguments["index"], File "/home/dray/anaconda3/envs/CI/lib/python3.10/site-packages/climate_indices/main.py", line 1324, in _parallel_process pool.map(_apply_along_axis_double, chunk_params) File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 364, in map return self._map_async(func, iterable, mapstar, chunksize).get() File "/home/dray/anaconda3/envs/CI/lib/python3.10/multiprocessing/pool.py", line 771, in get raise self._value ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments

On running: process_climate_indices --index spei --periodicity monthly --netcdf_precip LME_monthlyRAIN_AUSTNZ.nc --var_name_precip RAIN --netcdf_pet LME_monthlyEto_AUSTNZ.nc --var_name_pet PM_FAO_56 --output_file_base nclimgrid_SPEI12 --scales 12 --calibration_start_year 1986 --calibration_end_year 2005 --multiprocessing all

Link to rainfall and PET datasets: https://www.dropbox.com/scl/fi/8vpxhlyb7bllihlfm0th7/LME_monthlyRAIN_AUSTNZ.nc?dl=0&rlkey=wa2egqkal8oxpcmugyb46b647 https://www.dropbox.com/scl/fi/3kot9u44qpdt5bgp7v14p/LME_monthlyEto_AUSTNZ.nc?dl=0&rlkey=4uj8pgywhs1ffckuo3rrhqp21

Not sure if related or not but also get: NumbaDeprecationWarning: The 'nopython' keyword argument was not supplied to the 'numba.jit' decorator. The implicit default value for this argument is currently False, but it will be changed to True in Numba 0.59.0. See https://numba.readthedocs.io/en/stable/reference/deprecation.html#deprecation-of-object-mode-fall-back-behaviour-when-using-jit for details.

I have to run : export NUMBA_DISABLE_JIT=1 to get the job to even start to process to get to the L-moments error

I got the above command to run and worked a couple of weeks ago using version 1.10.12 but it wouldnt run this time so I updated to 1.10.13 and still not working.

conda list: _libgcc_mutex 0.1 conda_forge conda-forge _openmp_mutex 4.5 2_gnu conda-forge blosc 1.21.4 h0f2a231_0 conda-forge bzip2 1.0.8 h7f98852_4 conda-forge c-ares 1.19.1 hd590300_0 conda-forge ca-certificates 2023.5.7 hbcca054_0 conda-forge certifi 2023.5.7 pyhd8ed1ab_0 conda-forge cftime 1.6.2 py310hde88566_1 conda-forge click 8.1.3 pypi_0 pypi climate-indices 1.0.13 pypi_0 pypi cloudpickle 2.2.1 pypi_0 pypi curl 8.1.2 h409715c_0 conda-forge dask 2023.6.0 pypi_0 pypi esmf 8.4.2 nompi_h20110ff_0 conda-forge expat 2.5.0 hcb278e6_1 conda-forge fsspec 2023.6.0 pypi_0 pypi gsl 2.7 he838d99_0 conda-forge hdf4 4.2.15 h501b40f_6 conda-forge hdf5 1.14.0 nompi_hb72d44e_103 conda-forge icu 72.1 hcb278e6_0 conda-forge importlib-metadata 6.7.0 pypi_0 pypi keyutils 1.6.1 h166bdaf_0 conda-forge krb5 1.20.1 h81ceb04_0 conda-forge ld_impl_linux-64 2.40 h41732ed_0 conda-forge libaec 1.0.6 hcb278e6_1 conda-forge libblas 3.9.0 17_linux64_openblas conda-forge libcblas 3.9.0 17_linux64_openblas conda-forge libcurl 8.1.2 h409715c_0 conda-forge libedit 3.1.20191231 he28a2e2_2 conda-forge libev 4.33 h516909a_1 conda-forge libexpat 2.5.0 hcb278e6_1 conda-forge libffi 3.4.2 h7f98852_5 conda-forge libgcc-ng 13.1.0 he5830b7_0 conda-forge libgfortran-ng 13.1.0 h69a702a_0 conda-forge libgfortran5 13.1.0 h15d22d2_0 conda-forge libgomp 13.1.0 he5830b7_0 conda-forge libiconv 1.17 h166bdaf_0 conda-forge libjpeg-turbo 2.1.5.1 h0b41bf4_0 conda-forge liblapack 3.9.0 17_linux64_openblas conda-forge libnetcdf 4.9.2 nompi_h0f3d0bb_105 conda-forge libnghttp2 1.52.0 h61bc06f_0 conda-forge libnsl 2.0.0 h7f98852_0 conda-forge libopenblas 0.3.23 pthreads_h80387f5_0 conda-forge libsqlite 3.42.0 h2797004_0 conda-forge libssh2 1.11.0 h0841786_0 conda-forge libstdcxx-ng 13.1.0 hfd8a6a1_0 conda-forge libuuid 2.38.1 h0b41bf4_0 conda-forge libxml2 2.11.4 h0d562d8_0 conda-forge libzip 1.9.2 hc929e4a_1 conda-forge libzlib 1.2.13 hd590300_5 conda-forge llvmlite 0.40.1rc1 pypi_0 pypi locket 1.0.0 pypi_0 pypi lz4-c 1.9.4 hcb278e6_0 conda-forge nco 5.1.6 hd62b316_0 conda-forge ncurses 6.4 hcb278e6_0 conda-forge netcdf-fortran 4.6.1 nompi_h4f3791c_100 conda-forge netcdf4 1.6.4 nompi_py310hde23a83_100 conda-forge numba 0.57.0 pypi_0 pypi numpy 1.24.3 pypi_0 pypi openssl 3.1.1 hd590300_1 conda-forge packaging 23.1 pypi_0 pypi pandas 2.0.2 pypi_0 pypi partd 1.4.0 pypi_0 pypi pip 23.1.2 pyhd8ed1ab_0 conda-forge python 3.10.0 h543edf9_3_cpython conda-forge python-dateutil 2.8.2 pypi_0 pypi python_abi 3.10 3_cp310 conda-forge pytz 2023.3 pypi_0 pypi pyyaml 6.0 pypi_0 pypi readline 8.2 h8228510_1 conda-forge scipy 1.10.1 pypi_0 pypi setuptools 67.7.2 pyhd8ed1ab_0 conda-forge six 1.16.0 pypi_0 pypi snappy 1.1.10 h9fff704_0 conda-forge sqlite 3.42.0 h2c6b66d_0 conda-forge tempest-remap 2.2.0 h43474b4_0 conda-forge tk 8.6.12 h27826a3_0 conda-forge toolz 0.12.0 pypi_0 pypi tzdata 2023.3 pypi_0 pypi udunits2 2.2.28 hc3e0081_0 conda-forge wheel 0.40.0 pyhd8ed1ab_0 conda-forge xarray 2023.5.0 pypi_0 pypi xz 5.2.6 h166bdaf_0 conda-forge zipp 3.15.0 pypi_0 pypi zstd 1.5.2 h3eb15da_6 conda-forge

[monocongo]

There is an imminent update with code that I've tested against these Australian datasets, looking good now using the development branch named issue_522_pyproject_poetry in case you want to try it out right away.

[Darren-Ray]

Hi...thanks for that.

I did an uninstall and reinstall from the issue_522_pyproject_poetry tree.

The line 10 import typing in main.py needs to import Dict and Any as well.

When I fixed that, it ran, but again threw the L-moments error on Pearson typeIII :(

Is there a way to just run the Gamma distribution by the way?

And, from your answer, I assume the package is okay with cftime conventions? I have run the example .nc grids and they work fine by the way

[monocongo]

Thanks @Darren-Ray I appreciate the beta testing!

On my copy of that branch, I am having good results running the processing command. I've created a fresh conda environment and installed from source, like so:

conda create -n myvenv poetry pytest
conda activate myvenv
python -m poetry install
python -m poetry run pytest
cd ~/tmp/data
process_climate_indices --index spei --periodicity monthly --netcdf_precip LME_monthlyRAIN_AUSTNZ.nc --var_name_precip RAIN --netcdf_pet LME_monthlyEto_AUSTNZ.nc --var_name_pet PM_FAO_56 --output_file_base nclimgrid_SPEI12 --scales 12 --calibration_start_year 1986 --calibration_end_year 2005 --multiprocessing all

I have added cftime as a dependency which got me past an error I saw around that. So we should be OK for cftime conventions, but I'm not well-versed in that and can't promise too much there.

We used to have an option that allowed users to specify a single distribution but looking at the docs it seems that that was removed at some point, or maybe I am misremembering. Anyway, it shouldn't be too hard to re-implement that by adding a --distfit option or something like that.

[ben29med]

Hello every one, I tested the package "Another one" with fixing some problems and it's working fine. I would thank all of you, specially Mr. James Adams. All my good wishes to you.

Le sam. 1 juil. 2023 à 15:06, James Adams @.***> a écrit :

Thanks @Darren-Ray https://github.com/Darren-Ray I appreciate the beta testing!

On my copy of that branch, I am having good results running the processing command. I've created a fresh conda environment and installed from source, like so:

conda create -n myvenv poetry pytest conda activate myvenv python -m poetry install python -m poetry run pytest cd ~/tmp/data process_climate_indices --index spei --periodicity monthly --netcdf_precip LME_monthlyRAIN_AUSTNZ.nc --var_name_precip RAIN --netcdf_pet LME_monthlyEto_AUSTNZ.nc --var_name_pet PM_FAO_56 --output_file_base nclimgrid_SPEI12 --scales 12 --calibration_start_year 1986 --calibration_end_year 2005 --multiprocessing all

I have added cftime as a dependency which got me past an error I saw around that. So we should be OK for cftime conventions, but I'm not well-versed in that and can't promise too much there.

We used to have an option that allowed users to specify a single distribution but looking at the docs it seems that that was removed at some point, or maybe I am misremembering. Anyway, it shouldn't be too hard to re-implement that by adding a --distfit option or something like that.

— Reply to this email directly, view it on GitHub https://github.com/monocongo/climate_indices/issues/512#issuecomment-1615933063, or unsubscribe https://github.com/notifications/unsubscribe-auth/AJ5OUU6RAFMIARDVT2YKSPDXOAVEVANCNFSM6AAAAAAXVXSRYU . You are receiving this because you authored the thread.Message ID: @.***>

[chaowats]

Hi! @monocongo I have the same error, ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments, when calculating spi using data from CRU (https://crudata.uea.ac.uk/cru/data/hrg/cru_ts_4.07/cruts.2304141047.v4.07/pre/) I try to install it from the master branch and update Scipy but it is still not working. Any suggestion?

The example data from nclimgrid and cmorph is working well. The data from CHIRPS and ERA5 is working after reordering the variables and editing attributes.

[monocongo]

Thanks, @chaowats

Please confirm that you get the same result if installing this package into a fresh virtual environment using poetry. Please install from PyPI rather than github.

[ben29med]

@chaowats You can tell me how you installed the climate indices package because I'm traveling with the CHIRPS data. I installed the package and reorder the lat lon time, but it still not working !

[chaowats]

@monocongo Yes I did install to a new environment with the 2.0 version and master branch using pip install climate-indices from PyPI (but with conda environment, is that okay?) I got both the same Pearson Type III error.

[chaowats]

@ben29med I install using pip install climate-indices from PyPI but with the CHIRPS data I change precipitation units from mm/month to mm, longitude to lon, latitude to lat, and reorder to lat,lon,time before process_climate_indices --index spi.

[ben29med]

@chaowats please, send me the nclimgrid dataset example. Thank's !

[ZhiqiangD]

我在计算月尺度SPI Pearson 也碰到了同样的问题，可能是由于时间序列数据的问题，我将其加上很小的随机数就可计算。 nc_pre[yindex[i],xindex[i],:] = data + np.random.random(data.size)*0.01

[Naiver-008]

Hi! @monocongo I have the same error, ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments, when calculating spi using data from CRU (https://crudata.uea.ac.uk/cru/data/hrg/cru_ts_4.07/cruts.2304141047.v4.07/pre/) I try to install it from the master branch and update Scipy but it is still not working. Any suggestion?

The example data from nclimgrid and cmorph is working well. The data from CHIRPS and ERA5 is working after reordering the variables and editing attributes.

Hello, I came across the same problem with CRU data but ERA5 works fine. Did you find the solution for CRU data?

[rabons]

I also experience the L-moments error. I use the UCSB CHIRPS data. This follows a number of runtime warnings "divide by zero encountered in divide" and " invalid value encountered in multiply", which may be related to the L-moments error? Full log below. @monocongo is there a reference data set which should work?

The data can be found here: https://tubcloud.tu-berlin.de/s/4YcpA9RXLBXzka6

I am running the following commands:

 `$ conda create -n ci310 python=3.10
  $ conda activate ci310
  $ git clone https://github.com/monocongo/climate_indices.git
  $ cd climate_indices/
  $ poetry install
  $ process_climate_indices --index spi --periodicity monthly --netcdf_precip africa_chirps_1months_1981_2023.nc --var_name_precip precip --output_file_base CHIRPS --scales 3 --calibration_start_year 1981 --calibration_end_year 2023

Am I doing this right? Many thanks for you help, much appreciated. Below is the full terminal output.

2024-08-29  19:16:29 INFO Start time:    2024-08-29 19:16:29.147045
 2024-08-29  19:16:38 INFO Computing 3-month SPI/Pearson
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/miniconda3/envs/ci310/lib/python3.10/site-packages/scipy/stats/_distn_infrastructure.py:2176: RuntimeWarning: invalid value encountered in divide
  x = np.asarray((x - loc)/scale, dtype=dtyp)
2024-08-29  19:21:48 ERROR Unable to calculate Pearson Type III parameters due to invalid L-moments
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/miniconda3/envs/ci310/lib/python3.10/site-packages/scipy/stats/_distn_infrastructure.py:2176: RuntimeWarning: invalid value encountered in divide
  x = np.asarray((x - loc)/scale, dtype=dtyp)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:455: RuntimeWarning: divide by zero encountered in divide
  alpha = 4.0 / (skew * skew)
subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py:459: RuntimeWarning: invalid value encountered in multiply
  return loc - ((alpha * scale * skew) / 2.0)
subfolder/miniconda3/envs/ci310/lib/python3.10/site-packages/scipy/stats/_distn_infrastructure.py:2176: RuntimeWarning: invalid value encountered in divide
  x = np.asarray((x - loc)/scale, dtype=dtyp)
2024-08-29  19:31:17 ERROR Unable to calculate Pearson Type III parameters due to invalid L-moments
2024-08-29  19:31:17 ERROR Failed to complete
multiprocessing.pool.RemoteTraceback: 
"""
Traceback (most recent call last):
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 125, in worker
    result = (True, func(*args, **kwds))
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 48, in mapstar
    return list(map(*args))
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1292, in _apply_along_axis
    computed_array = np.apply_along_axis(func1d, axis=axis_index, arr=sub_array, parameters=args)
  File "<__array_function__ internals>", line 200, in apply_along_axis
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/site-packages/numpy/lib/shape_base.py", line 402, in apply_along_axis
    buff[ind] = asanyarray(func1d(inarr_view[ind], *args, **kwargs))
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1102, in _spi
    return indices.spi(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/indices.py", line 184, in spi
    values = compute.transform_fitted_pearson(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 624, in transform_fitted_pearson
    probabilities_of_zero, locs, scales, skews = pearson_parameters(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 307, in pearson_parameters
    prob, loc, scale, skew = calculate_time_step_params(time_step_values)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 256, in calculate_time_step_params
    params = lmoments.fit(time_step_values)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/lmoments.py", line 33, in fit
    raise ValueError(message)
ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments
"""

The above exception was the direct cause of the following exception:

Traceback (most recent call last):
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1617, in main
    _compute_write_index(kwrgs)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 997, in _compute_write_index
    _parallel_process(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1256, in _parallel_process
    pool.map(_apply_along_axis, chunk_params)
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 367, in map
    return self._map_async(func, iterable, mapstar, chunksize).get()
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 774, in get
    raise self._value
ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments
multiprocessing.pool.RemoteTraceback: 
"""
Traceback (most recent call last):
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 125, in worker
    result = (True, func(*args, **kwds))
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 48, in mapstar
    return list(map(*args))
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1292, in _apply_along_axis
    computed_array = np.apply_along_axis(func1d, axis=axis_index, arr=sub_array, parameters=args)
  File "<__array_function__ internals>", line 200, in apply_along_axis
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/site-packages/numpy/lib/shape_base.py", line 402, in apply_along_axis
    buff[ind] = asanyarray(func1d(inarr_view[ind], *args, **kwargs))
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1102, in _spi
    return indices.spi(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/indices.py", line 184, in spi
    values = compute.transform_fitted_pearson(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 624, in transform_fitted_pearson
    probabilities_of_zero, locs, scales, skews = pearson_parameters(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 307, in pearson_parameters
    prob, loc, scale, skew = calculate_time_step_params(time_step_values)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/compute.py", line 256, in calculate_time_step_params
    params = lmoments.fit(time_step_values)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/lmoments.py", line 33, in fit
    raise ValueError(message)
ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments
"""

The above exception was the direct cause of the following exception:

Traceback (most recent call last):
  File "subfolder/miniconda3/envs/ci310/bin/process_climate_indices", line 6, in <module>
    sys.exit(main())
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1617, in main
    _compute_write_index(kwrgs)
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 997, in _compute_write_index
    _parallel_process(
  File "subfolder/SPI/climate_indices/climate_indices/climate_indices/climate_indices/src/climate_indices/__main__.py", line 1256, in _parallel_process
    pool.map(_apply_along_axis, chunk_params)
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 367, in map
    return self._map_async(func, iterable, mapstar, chunksize).get()
  File "subfolder/miniconda3/envs/ci310/lib/python3.10/multiprocessing/pool.py", line 774, in get
    raise self._value
ValueError: Unable to calculate Pearson Type III parameters due to invalid L-moments`

[monocongo]

@rabons I see the same error when I run on my laptop. Looking at the data with ncdump -h shows that the time values are using a periodicity of "dekad":

        float precip(lat, lon, time) ;
                precip:_FillValue = -9999.f ;
                precip:units = "mm" ;
                precip:standard_name = "convective precipitation rate" ;
                precip:long_name = "Climate Hazards group InfraRed Precipitation with Stations" ;
                precip:time_step = "dekad" ;
                precip:missing_value = -9999.f ;
                precip:geospatial_lat_min = "NaN" ;
                precip:geospatial_lat_max = "NaN" ;
                precip:geospatial_lon_min = "NaN" ;
                precip:geospatial_lon_max = "NaN" ;
                precip:grid_mapping = "crs" ;

I was only recently made aware of this time step period but apparently, it's used by some researchers, so I'd like to support it as a new Periodicity enumerated type. I suspect the issue concerns a mismatch between monthly and dekad time steps for this error case. For example, we have two period (time step) types we support -- monthly and daily, and this will require adding a third type and modifying the processing loops to account for 36.5 dekads per year. We needed to do this before to add support for daily data (various modifications to allow for 366-day years) so this will be another iteration of that kind of refactoring.

[rabons]

@monocongo thank you for lookig into this!

time_step = "dekad" is actually a labeling error in this specific file, I think. Will double-check tomorrow. Would this label affect the SPI calculations in any way? I can re-produce it with a correct label.

EDIT: I have changed the file with correct labeling, it doesnt change the error with climate-indices. I have created this file with this script (slightly modified): https://bennyistanto.github.io/spi/chirpstif/

I actually do need a "dekadal" time step, but more important is that I make the monthly calculations work.

[monocongo]

OK then this isn't caused by the dekad time step since it's just a label. Thanks for the clarification @rabons

This error is raised in the lmoments.py code, here:

    # ensure the validity of the L-moments
    if (lmoments[1] <= 0) or (t3 >= 1):
        message = "Unable to calculate Pearson Type III parameters due to invalid L-moments"
        _logger.error(message)
        raise ValueError(message)

We can add some debugging around that to see if there are anomalies in the data etc. causing this to happen.

[rabons]

Thank you.

EDIT: SPI calculations seem to be working for a file clipped to a shapefile (in this case, Uganda boundaries). There may have been an unrecognized problem with my previous input file. The working input file: https://tubcloud.tu-berlin.de/s/omStM55dAmC2XAw and the spi result: https://tubcloud.tu-berlin.de/s/BE9E3ErefZD2MTn

The issue is hence solved for me when first clipping the input file (CHIRPS africa_monthly tiffs). Not sure what this is about. A few 'division by zero' errors persist, but I can still produce valid SPI output.