diku-dk / bfast Goto Github PK

I checked the file https://github.com/diku-dk/bfast/blob/master/bfast/monitor/python/base.py and the n_chunks parameter of the fit function is never used.
Do you prefer to implement a chunk loop or discard the parameter ?

Whenever one of the dates fed to model.fit() has the hour or minute set to something else then 0, the error
IndexError: boolean index did not match indexed array along dimension 1; dimension is 151 but corresponding boolean dimension is 152 is thrown.

This works also for the example code, when one value is modified like:

i = 34
di = dates_example[i]
dates_example[i] = datetime(di.year, di.month, di.day, 0, 34)

The error occurs in map_indices

A fix workaround for me is to use date instead of datetime objects.
I can also try to come up with a better fix, but that would need some time.

I have some image time series of length 730 time steps. I tried with multiple examples but I still get -2 corresponds to not enough data for a pixel for all pixels.

Could you please help me diagnose the issue?

I tried different values for freq and k in the BFASTMonitor with and without trend. I am working with python backend.

1. What are the key parameters that are required to be tuned?
2. Is there a verbose option?
3. What is n_chunks in model.fit call?

Hi,
I just installed the BFAST package.
But it does not work neither in my PC or in Google colab. Both present different errors.
Also both fail in the :
model.fit(data, dates, n_chunks=5, nan_value=-32768)
function. So the dataset and crop_data_dates() funtion works well.

Please, I will love some help, thanks.

I tried to use the bfast documentation and it seems that is it sno longer up to date.
To crop the data you suggested to use :

from bfast.utils import crop_data_dates
start_hist = datetime(2002, 1, 1)
start_monitor = datetime(2010, 1, 1)
end_monitor = datetime(2018, 1, 1)
data, dates = crop_data_dates(data_orig, dates, start_hist, end_monitor)
print("First date: {}".format(dates[0]))
print("Last date: {}".format(dates[-1]))
print("Shape of data array: {}".format(data.shape))

But utils is not part of bfast anymore. I created a PR to correct the documentation, feel free to accept it or change the __init__.pyof the bfast module to integrate utils in it.

While installing, it gives this warning in this colab notebook:
WARNING: Discarding https://files.pythonhosted.org/packages/d6/14/9418fc5dacae355be2ea1ec3d410d5d7b9c4d3d47c8eac1ef10dea20fd7f/bfast-0.7.tar.gz#sha256=649b8a6882db13634819b46ac056275683353591da77497d712236c595db2b79 (from https://pypi.org/simple/bfast/). Command errored out with exit status 1: python setup.py egg_info Check the logs for full command output.

So I checked the log, and it says that pyopencl is not found.

2021-08-30T16:09:55,905 Found link https://files.pythonhosted.org/packages/26/f1/b319368e760c8ba1017aadbd6de59cebf0b8cc2e93a1acad41db39e82a84/bfast-0.6.tar.gz#sha256=0126479167a744c5deec9288c1eaafc7ecb11471d4a32e8725ef2b285a1bf178 (from https://pypi.org/simple/bfast/), version: 0.6
2021-08-30T16:09:55,906 Found link https://files.pythonhosted.org/packages/d6/14/9418fc5dacae355be2ea1ec3d410d5d7b9c4d3d47c8eac1ef10dea20fd7f/bfast-0.7.tar.gz#sha256=649b8a6882db13634819b46ac056275683353591da77497d712236c595db2b79 (from https://pypi.org/simple/bfast/), version: 0.7
2021-08-30T16:09:55,906 Skipping link: not a file: https://pypi.org/simple/bfast/
2021-08-30T16:09:55,906 Given no hashes to check 27 links for project 'bfast': discarding no candidates
2021-08-30T16:09:55,915 Collecting bfast
2021-08-30T16:09:55,915 Created temporary directory: /tmp/pip-unpack-92hap_30
2021-08-30T16:09:55,941 Downloading bfast-0.7.tar.gz (306 kB)
2021-08-30T16:09:56,000 Added bfast from https://files.pythonhosted.org/packages/d6/14/9418fc5dacae355be2ea1ec3d410d5d7b9c4d3d47c8eac1ef10dea20fd7f/bfast-0.7.tar.gz#sha256=649b8a6882db13634819b46ac056275683353591da77497d712236c595db2b79 to build tracker '/tmp/pip-req-tracker-6oz_n0kw'
2021-08-30T16:09:56,001 Running setup.py (path:/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/setup.py) egg_info for package bfast
2021-08-30T16:09:56,001 Created temporary directory: /tmp/pip-pip-egg-info-km7t62k_
2021-08-30T16:09:56,001 Running command python setup.py egg_info
2021-08-30T16:09:56,645 Traceback (most recent call last):
2021-08-30T16:09:56,645 File "", line 1, in
2021-08-30T16:09:56,645 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/setup.py", line 16, in
2021-08-30T16:09:56,645 import bfast
2021-08-30T16:09:56,645 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/bfast/init.py", line 1, in
2021-08-30T16:09:56,645 from .models import BFASTMonitor
2021-08-30T16:09:56,645 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/bfast/models.py", line 1, in
2021-08-30T16:09:56,646 from bfast.monitor import BFASTMonitorPython
2021-08-30T16:09:56,646 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/bfast/monitor/init.py", line 1, in
2021-08-30T16:09:56,646 from .opencl import BFASTMonitorOpenCL
2021-08-30T16:09:56,646 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/bfast/monitor/opencl/init.py", line 1, in
2021-08-30T16:09:56,646 from .base import BFASTMonitorOpenCL
2021-08-30T16:09:56,646 File "/tmp/pip-install-xs33_vzd/bfast_d9f2a4c123a648ea8b440a89844c2869/bfast/monitor/opencl/base.py", line 13, in
2021-08-30T16:09:56,646 import pyopencl
2021-08-30T16:09:56,646 ModuleNotFoundError: No module named 'pyopencl'
2021-08-30T16:09:56,705 WARNING: Discarding https://files.pythonhosted.org/packages/d6/14/9418fc5dacae355be2ea1ec3d410d5d7b9c4d3d47c8eac1ef10dea20fd7f/bfast-0.7.tar.gz#sha256=649b8a6882db13634819b46ac056275683353591da77497d712236c595db2b79 (from https://pypi.org/simple/bfast/). Command errored out with exit status 1: python setup.py egg_info Check the logs for full command output.
2021-08-30T16:09:56,728 Created temporary directory: /tmp/pip-unpack-awqjgxhj
2021-08-30T16:09:56,800 Downloading bfast-0.6-py3-none-any.whl (811 kB)
2021-08-30T16:09:56,830 Added bfast from https://files.pythonhosted.org/packages/1b/93/81f97ddcb4a7dbfb08158ff5cf8abbc6b2cab07094b883cb0871513cd463/bfast-0.6-py3-none-any.whl#sha256=420c3c278cacdc7ec74ff1c593b29eeca8ed23e785cf78ea004b624c5c869e2f to build tracker '/tmp/pip-req-tracker-6oz_n0kw'
2021-08-30T16:09:56,830 Removed bfast from https://files.pythonhosted.org/packages/1b/93/81f97ddcb4a7dbfb08158ff5cf8abbc6b2cab07094b883cb0871513cd463/bfast-0.6-py3-none-any.whl#sha256=420c3c278cacdc7ec74ff1c593b29eeca8ed23e785cf78ea004b624c5c869e2f from build tracker '/tmp/pip-req-tracker-6oz_n0kw'

But running import pyopencl in colab gives no error. Moreover, version 0.6 still gets installed, and runs correctly using the opencl backend.

I'm interested in using this package for a project mapping deforestation with radar imagery. I've noticed that some open issues indicate some differences have existed between the results generated by this package and the original R version. Tt seems like some have been fixed and some issues are still yet to be resolved. If the maintainers have a moment, I'd like to get a sense of where the package stands and how it's been tested so I can decide if I should use it for my upcoming project. Thanks for releasing this cool project!

I'm trying to run the example on my Macbook Pro (no GPU) to test. I'm new to pyOpenCL but I read on the docs that it supports both CPU and GPU.

It is important to note that OpenCL is not restricted to GPUs. In fact, no special hardware is required to use OpenCL for computation–your existing CPU is enough. On Linux or macOS, type:

conda install pocl

to install a CPU-based OpenCL driver. On Windows, you may install e.g. the CPU OpenCL driver from Intel. On macOS, pocl can offer a marked robustness (and, sometimes, performance) improvement over the OpenCL drivers built into the operating system.

I installed pyopencl from conda-forge:

pyopencl 2018.2.5 py37h9888f84_0 conda-forge

And first tried to run the example with this installation:

import os
import wget
import numpy
from datetime import datetime

# download and parse input data
ifile_meta = "../data/bfast/peru_small/dates.txt"
ifile_data = "../data/bfast/peru_small/data.npy"

if not os.path.isdir("../data/bfast/peru_small"):
    os.makedirs("../data/bfast/peru_small")

if not os.path.exists(ifile_meta):
    url = 'https://sid.erda.dk/share_redirect/fcwjD77gUY/dates.txt'
    wget.download(url, ifile_meta)
if not os.path.exists(ifile_data):
    url = 'https://sid.erda.dk/share_redirect/fcwjD77gUY/data.npy'
    wget.download(url, ifile_data)

data_orig = numpy.load(ifile_data)
with open(ifile_meta) as f:
    dates = f.read().split('\n')
    dates = [datetime.strptime(d, '%Y-%m-%d') for d in dates if len(d) > 0]


from bfast.utils import crop_data_dates
start_hist = datetime(2002, 1, 1)
start_monitor = datetime(2010, 1, 1)
end_monitor = datetime(2018, 1, 1)
data, dates = crop_data_dates(data_orig, dates, start_hist, end_monitor)
print("First date: {}".format(dates[0]))
print("Last date: {}".format(dates[-1]))
print("Shape of data array: {}".format(data.shape))

from bfast import BFASTMonitor

model = BFASTMonitor(
            start_monitor,
            freq=365,
            k=3,
            hfrac=0.25,
            trend=False,
            level=0.05,
            backend='opencl',
            device_id=0,
        )
model.fit(data, dates, n_chunks=5, nan_value=-32768)

print("Detected breaks")
# -2 corresponds to not enough data for a pixel
# -1 corresponds to "no breaks detected"
# idx with isx>=0 corresponds to the position of the first break
print(model.breaks)

I get this error

---------------------------------------------------------------------------
RuntimeError                              Traceback (most recent call last)
<ipython-input-1-6b81a4d59de0> in <module>
     45             device_id=0,
     46         )
---> 47 model.fit(data, dates, n_chunks=5, nan_value=-32768)
     48 
     49 print("Detected breaks")

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/bfast/models.py in fit(self, data, dates, n_chunks, nan_value)
    181                  verbose=self.verbose,
    182                  platform_id=self.platform_id,
--> 183                  device_id=self.device_id
    184                 )
    185 

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/bfast/monitor/opencl/base.py in __init__(self, start_monitor, freq, k, hfrac, trend, level, period, detailed_results, find_magnitudes, verbose, platform_id, device_id)
    110                                  default_tile_size=8,
    111                                  default_reg_tile_size=3,
--> 112                                  sizes=self._get_futhark_params())
    113 
    114     def _init_device(self, platform_id, device_id):

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/bfast/monitor/opencl/bfastfinal.py in __init__(self, command_queue, interactive, platform_pref, device_pref, default_group_size, default_num_groups, default_tile_size, default_reg_tile_size, default_threshold, sizes)
  48814                                                                            "value": None},
  48815                                         "mainMagnitude.tile_size_41656": {"class": "tile_size", "value": None},
> 48816                                         "remove_nans.segmap_group_size_29490": {"class": "group_size", "value": None}})
  48817     self.builtinzhreplicate_f32zireplicate_44709_var = program.builtinzhreplicate_f32zireplicate_44709
  48818     self.builtinzhreplicate_i32zireplicate_44718_var = program.builtinzhreplicate_i32zireplicate_44718

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/bfast/monitor/opencl/bfastfinal.py in initialise_opencl_object(self, program_src, command_queue, interactive, platform_pref, device_pref, default_group_size, default_num_groups, default_tile_size, default_reg_tile_size, default_threshold, size_heuristics, required_types, all_sizes, user_sizes)
    220         return cl.Program(self.ctx, program_src).build(
    221             ["-DLOCKSTEP_WIDTH={}".format(lockstep_width)]
--> 222             + ["-D{}={}".format(s.replace('z', 'zz').replace('.', 'zi').replace('#', 'zh'),v) for (s,v) in self.sizes.items()])
    223 
    224 def opencl_alloc(self, min_size, tag):

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/pyopencl/__init__.py in build(self, options, devices, cache_dir)
    508                         self._context, self._source, options_bytes, devices,
    509                         cache_dir=cache_dir, include_path=include_path),
--> 510                     options_bytes=options_bytes, source=self._source)
    511 
    512             if was_cached:

~/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/pyopencl/__init__.py in _build_and_catch_errors(self, build_func, options_bytes, source)
    552         # Python 3.2 outputs the whole list of currently active exceptions
    553         # This serves to remove one (redundant) level from that nesting.
--> 554         raise err
    555 
    556     # }}}

RuntimeError: clBuildProgram failed: BUILD_PROGRAM_FAILURE - clBuildProgram failed: BUILD_PROGRAM_FAILURE - clBuildProgram failed: BUILD_PROGRAM_FAILURE

Build on <pyopencl.Device 'Intel(R) Core(TM) i7-1068NG7 CPU @ 2.30GHz' on 'Apple' at 0x7fced9aea050>:


(options: -DLOCKSTEP_WIDTH=1 -Dbuiltinzhreplicate_f32zigroup_sizze_44712=32 -Dbuiltinzhreplicate_i32zigroup_sizze_44721=32 -DmainziRx_41192=3 -DmainziRx_41933=3 -DmainziRy_41193=3 -DmainziRy_41934=3 -DmainziTk_41189=8 -DmainziTk_41930=8 -DmainziTx_41044=8 -DmainziTx_41190=8 -DmainziTx_41931=8 -DmainziTy_41045=8 -DmainziTy_41191=8 -DmainziTy_41932=8 -Dmainzigroup_sizze_44243=32 -Dmainzisegmap_group_sizze_37619=32 -Dmainzisegmap_group_sizze_37797=32 -Dmainzisegmap_group_sizze_37925=32 -Dmainzisegmap_group_sizze_37957=32 -Dmainzisegmap_group_sizze_38004=32 -Dmainzisegmap_group_sizze_38489=32 -Dmainzisegmap_group_sizze_38654=32 -Dmainzisegmap_group_sizze_38708=32 -Dmainzisegmap_group_sizze_38775=32 -Dmainzisegmap_group_sizze_38869=32 -Dmainzisegmap_group_sizze_39049=32 -Dmainzisegmap_group_sizze_39190=32 -Dmainzisegmap_group_sizze_39322=32 -Dmainzisegmap_group_sizze_39603=32 -Dmainzisegmap_group_sizze_39678=32 -Dmainzisegmap_group_sizze_39827=32 -Dmainzisegmap_group_sizze_39929=32 -Dmainzisegmap_group_sizze_40076=32 -Dmainzisegmap_group_sizze_40200=32 -Dmainzisegmap_group_sizze_40570=32 -Dmainzisegmap_group_sizze_40712=32 -Dmainzisegmap_num_groups_37959=8 -Dmainzisegmap_num_groups_38006=8 -Dmainzisegmap_num_groups_39051=8 -Dmainzisegmap_num_groups_39192=8 -Dmainzisegmap_num_groups_39324=8 -Dmainzisegmap_num_groups_40714=8 -Dmainzisegred_group_sizze_38064=32 -Dmainzisegred_group_sizze_39111=32 -Dmainzisegred_group_sizze_39248=32 -Dmainzisegred_group_sizze_39378=32 -Dmainzisegred_group_sizze_39944=32 -Dmainzisegred_group_sizze_39965=32 -Dmainzisegred_group_sizze_40032=32 -Dmainzisegred_group_sizze_40116=32 -Dmainzisegred_group_sizze_40617=32 -Dmainzisegred_num_groups_38066=8 -Dmainzisegred_num_groups_39113=8 -Dmainzisegred_num_groups_39250=8 -Dmainzisegred_num_groups_39380=8 -Dmainzisegred_num_groups_39946=8 -Dmainzisegred_num_groups_39967=8 -Dmainzisegred_num_groups_40034=8 -Dmainzisegred_num_groups_40118=8 -Dmainzisegred_num_groups_40619=8 -Dmainzisegscan_group_sizze_39687=32 -Dmainzisegscan_group_sizze_40671=32 -Dmainzisegscan_num_groups_39689=8 -Dmainzisegscan_num_groups_40673=8 -Dmainzisuff_intra_par_11=128 -Dmainzisuff_intra_par_13=128 -Dmainzisuff_intra_par_24=235 -Dmainzisuff_intra_par_29=113 -Dmainzisuff_intra_par_34=122 -Dmainzisuff_outer_par_16=2000000000 -Dmainzisuff_outer_par_17=892448 -Dmainzisuff_outer_par_18=2000000000 -Dmainzisuff_outer_par_19=892448 -Dmainzisuff_outer_par_20=2000000000 -Dmainzisuff_outer_par_21=26215660 -Dmainzisuff_outer_par_28=2000000000 -Dmainzisuff_outer_par_31=111556 -Dmainzisuff_outer_par_6=2000000000 -Dmainzisuff_outer_par_7=2000000000 -Dmainzisuff_outer_par_8=7139584 -Dmainzitile_sizze_41656=8 -DmainDetailedziRx_41192=3 -DmainDetailedziRx_41933=3 -DmainDetailedziRy_41193=3 -DmainDetailedziRy_41934=3 -DmainDetailedziTk_41189=8 -DmainDetailedziTk_41930=8 -DmainDetailedziTx_41044=8 -DmainDetailedziTx_41190=8 -DmainDetailedziTx_41931=8 -DmainDetailedziTy_41045=8 -DmainDetailedziTy_41191=8 -DmainDetailedziTy_41932=8 -DmainDetailedzigroup_sizze_44256=32 -DmainDetailedzisegmap_group_sizze_29648=32 -DmainDetailedzisegmap_group_sizze_29826=32 -DmainDetailedzisegmap_group_sizze_29954=32 -DmainDetailedzisegmap_group_sizze_29986=32 -DmainDetailedzisegmap_group_sizze_30033=32 -DmainDetailedzisegmap_group_sizze_30518=32 -DmainDetailedzisegmap_group_sizze_30683=32 -DmainDetailedzisegmap_group_sizze_30737=32 -DmainDetailedzisegmap_group_sizze_30804=32 -DmainDetailedzisegmap_group_sizze_30898=32 -DmainDetailedzisegmap_group_sizze_31078=32 -DmainDetailedzisegmap_group_sizze_31219=32 -DmainDetailedzisegmap_group_sizze_31351=32 -DmainDetailedzisegmap_group_sizze_31632=32 -DmainDetailedzisegmap_group_sizze_31707=32 -DmainDetailedzisegmap_group_sizze_31856=32 -DmainDetailedzisegmap_group_sizze_31958=32 -DmainDetailedzisegmap_group_sizze_32105=32 -DmainDetailedzisegmap_group_sizze_32229=32 -DmainDetailedzisegmap_group_sizze_32480=32 -DmainDetailedzisegmap_group_sizze_32602=32 -DmainDetailedzisegmap_group_sizze_32659=32 -DmainDetailedzisegmap_group_sizze_33204=32 -DmainDetailedzisegmap_group_sizze_33256=32 -DmainDetailedzisegmap_group_sizze_33291=32 -DmainDetailedzisegmap_group_sizze_33412=32 -DmainDetailedzisegmap_num_groups_29988=8 -DmainDetailedzisegmap_num_groups_30035=8 -DmainDetailedzisegmap_num_groups_31080=8 -DmainDetailedzisegmap_num_groups_31221=8 -DmainDetailedzisegmap_num_groups_31353=8 -DmainDetailedzisegmap_num_groups_33414=8 -DmainDetailedzisegred_group_sizze_30093=32 -DmainDetailedzisegred_group_sizze_31140=32 -DmainDetailedzisegred_group_sizze_31277=32 -DmainDetailedzisegred_group_sizze_31407=32 -DmainDetailedzisegred_group_sizze_31973=32 -DmainDetailedzisegred_group_sizze_31994=32 -DmainDetailedzisegred_group_sizze_32061=32 -DmainDetailedzisegred_group_sizze_32145=32 -DmainDetailedzisegred_group_sizze_33317=32 -DmainDetailedzisegred_num_groups_30095=8 -DmainDetailedzisegred_num_groups_31142=8 -DmainDetailedzisegred_num_groups_31279=8 -DmainDetailedzisegred_num_groups_31409=8 -DmainDetailedzisegred_num_groups_31975=8 -DmainDetailedzisegred_num_groups_31996=8 -DmainDetailedzisegred_num_groups_32063=8 -DmainDetailedzisegred_num_groups_32147=8 -DmainDetailedzisegred_num_groups_33319=8 -DmainDetailedzisegscan_group_sizze_31716=32 -DmainDetailedzisegscan_group_sizze_33371=32 -DmainDetailedzisegscan_num_groups_31718=8 -DmainDetailedzisegscan_num_groups_33373=8 -DmainDetailedzisuff_intra_par_11=32 -DmainDetailedzisuff_intra_par_13=32 -DmainDetailedzisuff_intra_par_24=32 -DmainDetailedzisuff_intra_par_29=32 -DmainDetailedzisuff_intra_par_37=32 -DmainDetailedzisuff_outer_par_16=8 -DmainDetailedzisuff_outer_par_17=8 -DmainDetailedzisuff_outer_par_18=8 -DmainDetailedzisuff_outer_par_19=8 -DmainDetailedzisuff_outer_par_20=8 -DmainDetailedzisuff_outer_par_21=8 -DmainDetailedzisuff_outer_par_28=8 -DmainDetailedzisuff_outer_par_31=8 -DmainDetailedzisuff_outer_par_6=8 -DmainDetailedzisuff_outer_par_7=8 -DmainDetailedzisuff_outer_par_8=8 -DmainDetailedzitile_sizze_41656=8 -DmainMagnitudeziRx_41192=3 -DmainMagnitudeziRx_41933=3 -DmainMagnitudeziRy_41193=3 -DmainMagnitudeziRy_41934=3 -DmainMagnitudeziTk_41189=8 -DmainMagnitudeziTk_41930=8 -DmainMagnitudeziTx_41044=8 -DmainMagnitudeziTx_41190=8 -DmainMagnitudeziTx_41931=8 -DmainMagnitudeziTy_41045=8 -DmainMagnitudeziTy_41191=8 -DmainMagnitudeziTy_41932=8 -DmainMagnitudezigroup_sizze_44244=32 -DmainMagnitudezisegmap_group_sizze_33713=32 -DmainMagnitudezisegmap_group_sizze_33891=32 -DmainMagnitudezisegmap_group_sizze_34019=32 -DmainMagnitudezisegmap_group_sizze_34051=32 -DmainMagnitudezisegmap_group_sizze_34098=32 -DmainMagnitudezisegmap_group_sizze_34583=32 -DmainMagnitudezisegmap_group_sizze_34748=32 -DmainMagnitudezisegmap_group_sizze_34802=32 -DmainMagnitudezisegmap_group_sizze_34869=32 -DmainMagnitudezisegmap_group_sizze_34963=32 -DmainMagnitudezisegmap_group_sizze_35143=32 -DmainMagnitudezisegmap_group_sizze_35284=32 -DmainMagnitudezisegmap_group_sizze_35416=32 -DmainMagnitudezisegmap_group_sizze_35697=32 -DmainMagnitudezisegmap_group_sizze_35772=32 -DmainMagnitudezisegmap_group_sizze_35921=32 -DmainMagnitudezisegmap_group_sizze_36023=32 -DmainMagnitudezisegmap_group_sizze_36170=32 -DmainMagnitudezisegmap_group_sizze_36294=32 -DmainMagnitudezisegmap_group_sizze_36545=32 -DmainMagnitudezisegmap_group_sizze_36667=32 -DmainMagnitudezisegmap_group_sizze_36724=32 -DmainMagnitudezisegmap_group_sizze_37222=32 -DmainMagnitudezisegmap_group_sizze_37364=32 -DmainMagnitudezisegmap_num_groups_34053=8 -DmainMagnitudezisegmap_num_groups_34100=8 -DmainMagnitudezisegmap_num_groups_35145=8 -DmainMagnitudezisegmap_num_groups_35286=8 -DmainMagnitudezisegmap_num_groups_35418=8 -DmainMagnitudezisegmap_num_groups_37366=8 -DmainMagnitudezisegred_group_sizze_34158=32 -DmainMagnitudezisegred_group_sizze_35205=32 -DmainMagnitudezisegred_group_sizze_35342=32 -DmainMagnitudezisegred_group_sizze_35472=32 -DmainMagnitudezisegred_group_sizze_36038=32 -DmainMagnitudezisegred_group_sizze_36059=32 -DmainMagnitudezisegred_group_sizze_36126=32 -DmainMagnitudezisegred_group_sizze_36210=32 -DmainMagnitudezisegred_group_sizze_37269=32 -DmainMagnitudezisegred_num_groups_34160=8 -DmainMagnitudezisegred_num_groups_35207=8 -DmainMagnitudezisegred_num_groups_35344=8 -DmainMagnitudezisegred_num_groups_35474=8 -DmainMagnitudezisegred_num_groups_36040=8 -DmainMagnitudezisegred_num_groups_36061=8 -DmainMagnitudezisegred_num_groups_36128=8 -DmainMagnitudezisegred_num_groups_36212=8 -DmainMagnitudezisegred_num_groups_37271=8 -DmainMagnitudezisegscan_group_sizze_35781=32 -DmainMagnitudezisegscan_group_sizze_37323=32 -DmainMagnitudezisegscan_num_groups_35783=8 -DmainMagnitudezisegscan_num_groups_37325=8 -DmainMagnitudezisuff_intra_par_11=128 -DmainMagnitudezisuff_intra_par_13=2000000000 -DmainMagnitudezisuff_intra_par_24=235 -DmainMagnitudezisuff_intra_par_29=113 -DmainMagnitudezisuff_intra_par_37=122 -DmainMagnitudezisuff_outer_par_16=2000000000 -DmainMagnitudezisuff_outer_par_17=892448 -DmainMagnitudezisuff_outer_par_18=111556 -DmainMagnitudezisuff_outer_par_19=2000000000 -DmainMagnitudezisuff_outer_par_20=2000000000 -DmainMagnitudezisuff_outer_par_21=26215660 -DmainMagnitudezisuff_outer_par_28=2000000000 -DmainMagnitudezisuff_outer_par_31=111556 -DmainMagnitudezisuff_outer_par_6=2000000000 -DmainMagnitudezisuff_outer_par_7=2000000000 -DmainMagnitudezisuff_outer_par_8=7139584 -DmainMagnitudezitile_sizze_41656=8 -Dremove_nanszisegmap_group_sizze_29490=32 -I /Users/rave/opt/miniconda3/envs/pybayts/lib/python3.7/site-packages/pyopencl/cl)
(source saved as /var/folders/t5/9x0zlgbs2dz25fcjxvh9fpym0000gn/T/tmpqmd1y_4f.cl)

I then saw conda-forge/pyopencl-feedstock#26 (comment) so I then installed pocl in my conda environment with

conda install -c conda-forge osx-pocl-opencl pocl pyopencl==2018.2.5

But when I ran the example again, I got a dead kernel error (I'm running it in a jupyter notebook).

Any tips on how to set up the environment on a CPU-only machine?

Hello,

I went throw the documentation and it seems that we can only get the date of the first break for one pixel. Is it possible to get the full list of breaks (in the case where there are several) for one pixel during the monitoring period ?

Thank you very much

Can only the first break be discovered using this package or does it also support the detection of multiple breaks in the time series? If yes, how could it be done, as the example code (at least as I understand) only detects the first break.

Can we amortize the data transfer overheads (host -> device and device -> host with kernel execution) when running with multiple chunks using multithreading?

Hello Dmitry!

Thanks for the new commits! 👍

In the meantime, I found something that I think you know, but I would like to share.

For pixels that have few useful time samples, as you know, in Python they are treated with the condition:

From some results I obtained with Tile_1 data, I found some discrepancies in the results for magnitudes in Python vs OpenCL:

Where magnsPB is magnitudes in Python and magnsCLB is magnitudes in OpenCL, with a MSE (mean squared error) of 142647860.0.

This is the plot for the magnitudes <=0:

As we know, the algorithms are different so the results, but I managed to get some better results when I applied the same condition in OpenCL, after returning detailed results. Something like:

Namely, for this particular example, I like to use threshold 9 instead of 5. I obtained the results:

With a MSE of 0.30474687.

Do you think it is a good idea for the user to choose the value for the condition?

It's just an idea, not a big problem if we have good data to process.

Thank you!

Thanks for your work. Are there plans to integrate the "traditional" bfast() function as well?

Hi,

Thanks for the package, I am just wondering could I apply the algorithm on 1-d time series data? I checked the format of input data which should be (N,W,H) only for image data. Thank you!

issue

I wanted to use the specific develop branch in my requirement file as it's currently the most advanced one but it cannot be installed this way.

Reproductible example:

create a requirements.txt file containing the following:

git+git://github.com/diku-dk/bfast.git@develop#egg=bfast

then create a fresh venv anywhere in your home directory and activate it:

python3 -m venv test
source test/bin/activate 

and try to recursively install the previously created reqs:

pip3 install -r requirements.txt

You will get the following error:

Collecting bfast from git+git://github.com/diku-dk/bfast.git@develop#egg=bfast (from -r ../license.txt (line 1))
Cloning git://github.com/diku-dk/bfast.git (to develop) to /tmp/pip-build-7v8ihvfi/bfast
Complete output from command python setup.py egg_info:
Traceback (most recent call last):
File "", line 1, in
File "/tmp/pip-build-7v8ihvfi/bfast/setup.py", line 16, in
import bfast
File "/tmp/pip-build-7v8ihvfi/bfast/bfast/init.py", line 1, in
from .models import BFASTMonitor
File "/tmp/pip-build-7v8ihvfi/bfast/bfast/models.py", line 1, in
from bfast.monitor import BFASTMonitorPython
File "/tmp/pip-build-7v8ihvfi/bfast/bfast/monitor/init.py", line 1, in
from .opencl import BFASTMonitorOpenCL
File "/tmp/pip-build-7v8ihvfi/bfast/bfast/monitor/opencl/init.py", line 1, in
from .base import BFASTMonitorOpenCL
File "/tmp/pip-build-7v8ihvfi/bfast/bfast/monitor/opencl/base.py", line 11, in
import pandas
ModuleNotFoundError: No module named 'pandas'

---

Command "python setup.py egg_info" failed with error code 1 in /tmp/pip-build-7v8ihvfi/bfast/

problem

The issue comes from the __init__.py file of bfast.
In setup.py you import bfast to retrieve the version number. This action is not possible if the other libs are not yet installed because there is an import of bfastMonitor in the __init__.py file.It's the case if you use a -r installation process (pip install them all at once and not one after another).

I would suggest to code the version number in a separated file if you really want to import it or to hard code it in the setup.py file

Hello, do you have any idea when the magnitude calculation for the openCL version will be available?

I'm using the y_error and the ns values that are returned from bfastfinaldetailed, but it takes some time to do the calculation on the CPU, as well as a large allocation of GPU memory necessary to perform the bfastfinaldetailed. I have already looked at the bfastfinal and bfastfinaldetailed files a lot, trying to understand which kernels do what, but they are very unreadable.
I saw that you made commits to the opencl-backend-debug branch with work for magnitude calculation for the openCL, but I was unable to use the updates.

Thank you!

Hi,

I had a look at the python version and if I'm not mistaken a stable history period is never determined as in section 2.4 of the paper.

Am I correct, or did I misunderstand the code? I tried to read the futhark implementation and it doesn't seem to be there either.

Hey there !

Thanks for the python implementation of the bfast algo !
However, I'm still a bit confused about the results of the "model.breaks" array.
What are the meanings of the values "-2", "-1", etc.
I guess that the positive values corresponds to the index of the image where the break is detected ?

Thanks for the clarifications !

Hello,

It's great we now have bfast for python.
I'm quite new in the field. I'm trying to implement the technic in a time-series with 180 dates, but I'm getting the value error that the length of the passed values is not equal to the required, after processing row 0. I've also tested the example "small-Peru" and raised the same error. I've checked that the required dependencies meet the desired versions and also reinstalled NumPy, but I keep getting the same error. Below are the errors I get.

Small Peru Example

`First date: 2003-07-14 00:00:00
Last date: 2016-09-27 00:00:00
Shape of data array: (216, 334, 334)
Processing row 0
Traceback (most recent call last):
File "C:/Users/vfota/diploma/test_bfast.py", line 63, in
model.fit(data, dates, nan_value=-32768)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\models.py", line 180, in fit
nan_value=nan_value)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 144, in fit
self.fit_single(y, dates)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 175, in fit_single
self.mapped_indices = self._map_indices(dates).astype(numpy.int32)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 299, in _map_indices
indices = numpy.argwhere(~numpy.isnan(ts)).T[0]
File "<array_function internals>", line 6, in argwhere
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\numeric.py", line 584, in argwhere
return transpose(nonzero(a))
File "<array_function internals>", line 6, in nonzero
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 1896, in nonzero
return _wrapfunc(a, 'nonzero')
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 58, in _wrapfunc
return _wrapit(obj, method, *args, **kwds)
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 51, in _wrapit
result = wrap(result)
File "C:\Users\vfota\diploma\venv\lib\site-packages\pandas\core\generic.py", line 1917, in array_wrap
return self._constructor(result, **d).finalize(self)
File "C:\Users\vfota\diploma\venv\lib\site-packages\pandas\core\series.py", line 292, in init
f"Length of passed values is {len(data)}, "
ValueError: Length of passed values is 1, index implies 5114.

Process finished with exit code 1
`

My Time-Series

`First date: 2000-01-21 00:00:00
Last date: 2015-12-24 00:00:00
Shape of data array: (181, 2669, 4094)
Start monitor: 2010-02-17 00:00:00
Processing row 0
Traceback (most recent call last):
File "C:/Users/vfota/diploma/bfast_euboea.py", line 71, in
nan_value=np.nan)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\models.py", line 180, in fit
nan_value=nan_value)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 144, in fit
self.fit_single(y, dates)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 175, in fit_single
self.mapped_indices = self._map_indices(dates).astype(numpy.int32)
File "C:\Users\vfota\diploma\venv\lib\site-packages\bfast\monitor\python\base.py", line 299, in _map_indices
indices = numpy.argwhere(~numpy.isnan(ts)).T[0]
File "<array_function internals>", line 6, in argwhere
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\numeric.py", line 584, in argwhere
return transpose(nonzero(a))
File "<array_function internals>", line 6, in nonzero
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 1896, in nonzero
return _wrapfunc(a, 'nonzero')
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 58, in _wrapfunc
return _wrapit(obj, method, *args, **kwds)
File "C:\Users\vfota\diploma\venv\lib\site-packages\numpy\core\fromnumeric.py", line 51, in _wrapit
result = wrap(result)
File "C:\Users\vfota\diploma\venv\lib\site-packages\pandas\core\generic.py", line 1917, in array_wrap
return self._constructor(result, **d).finalize(self)
File "C:\Users\vfota\diploma\venv\lib\site-packages\pandas\core\series.py", line 292, in init
f"Length of passed values is {len(data)}, "
ValueError: Length of passed values is 1, index implies 5844.

Process finished with exit code 1`

bfast==0.5
matplotlib==3.1.3
numpy==1.19.2
pandas==1.0.0
pyopencl== 2020.2.2+cl12
pyOpenSSL== 19.0.0
scikit-learn==0.23.2
scipy==1.5.2
Sphinx==3.2.1
sphinx-bootstrap-theme==0.7.1
wget== 3.2
numpydoc==1.0.0

Good morning, it is possible to obtain different breaks and means values for the same data (and dates) inputs, when I run the algorithm in python and openCL?

The input data is a simple one with 124 series 17x17. In this case, I executed the algorithm in python, and next, I executed in openCL and obtained different values for breaks and means. I used in both cases the default values for k = 3, freq = 365, trend = False, hfrac = 0.25, level = 0.05 and I just changed the backend value. I must say that for some start_monitor dates the results for the breaks are equal (which is good), but not always. It is possible to obtain different breaks and means values when I change only the backend value, or I'm doing something wrong?

Thank you all, in advance!
Best regards

There are problems with tuning thresholds in base.py which make it perform much slower (up to 5x speedup is possible with proper tuning). Need to look into retuning the code with the newest version of the Futhark compiler

I would create one instance of BFASTMonitor for each GPU (you can use device_id to choose GPU), then run fit() on each window. Alternatively, you can use bigger windows, and use n_chunks, to let BFASTMonitor split the window into smaller chunks. Depends on the size of the datasets that you are using and the available RAM.

I checked the BFASTMonitor object and _model is not set once and for all when BFASTMonitor is instantiated but every time the fit method is called. Meaning that for GPU implementation the setting of the device is reset every time.

I don't see any parameter relative to the _model in the the fit method, wouldn't it be more semantic to move it to the __init__ method ?

Happy to create a PR

I am using the following Python and R bfast versions at the date of the creation of this issue.

git clone --single-branch --branch develop https://github.com/gieseke/bfast.git
git clone --single-branch --branch master https://github.com/bfast2/bfast.git

To minimize complexity, I created these 2 self-contained python and R scripts, that are almost direct translations:

import datetime
import numpy as np
import bisect
import bfast

indexedRasterTimeSeries = np.array([[[3141]],[[0]],[[2153]],[[0]],[[2791]],[[2795]],[[2109]],[[2139]],[[2077]],[[2200]],[[0]],[[1260]],[[1301]],[[2496]],[[1429]],[[0]],[[0]],[[0]],[[0]],[[3122]],[[2446]],[[2428]],[[1647]],[[1266]],[[0]],[[25]],[[-996]],[[-21]],[[0]],[[0]],[[0]],[[0]],[[466]],[[1454]],[[0]],[[0]],[[2084]],[[2600]],[[0]],[[2133]],[[1865]],[[0]],[[596]],[[410]],[[0]],[[0]],[[0]],[[0]],[[2438]],[[0]],[[0]],[[0]],[[0]],[[3382]],[[3753]],[[4647]],[[0]],[[2989]],[[3356]],[[3223]],[[2732]],[[3536]],[[2332]],[[0]],[[1957]],[[1325]],[[0]],[[1990]],[[0]],[[0]],[[2740]],[[3510]],[[0]],[[0]],[[3019]],[[0]],[[3676]],[[3790]],[[3370]],[[0]],[[3504]],[[0]],[[0]],[[2984]],[[2762]],[[2794]],[[2845]],[[0]],[[0]],[[3076]],[[0]],[[3213]],[[0]],[[0]],[[0]],[[0]],[[0]],[[4071]],[[0]],[[0]],[[0]],[[0]],[[3387]],[[0]],[[3456]],[[3222]],[[3173]],[[2609]],[[2650]],[[2535]],[[2441]],[[2191]],[[2472]],[[2350]],[[2847]],[[0]],[[0]],[[3484]],[[0]],[[4704]],[[3412]],[[3687]],[[0]],[[0]],[[0]],[[0]],[[3946]],[[0]],[[0]],[[0]],[[0]],[[3764]],[[3266]],[[3451]],[[3552]],[[3289]],[[2517]],[[2588]],[[0]],[[0]],[[1487]],[[1746]],[[2084]],[[1884]],[[1917]],[[0]],[[3765]],[[0]],[[0]],[[0]],[[0]],[[3950]],[[3663]],[[3733]],[[4007]],[[4202]],[[3566]],[[3926]],[[3573]],[[3804]],[[3504]],[[3678]],[[3318]],[[3255]],[[0]],[[0]],[[2768]],[[3016]],[[0]],[[0]],[[3641]],[[3685]],[[4161]],[[0]],[[3975]],[[0]],[[0]],[[4019]],[[0]],[[4026]],[[3697]],[[4053]],[[4292]],[[3741]],[[3911]],[[3575]],[[3467]],[[3410]],[[0]],[[3291]],[[3473]],[[0]],[[0]],[[3467]],[[0]],[[3801]],[[3022]],[[0]],[[0]],[[3363]],[[3742]],[[4035]],[[0]],[[0]],[[0]],[[3996]],[[0]],[[0]],[[3848]],[[0]],[[0]],[[3709]],[[3312]],[[3148]],[[2734]],[[2944]],[[0]],[[0]],[[0]],[[0]],[[4086]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[0]],[[4200]]])
observationDates = ['2010-02-13', '2010-03-09', '2010-04-10', '2010-04-18', '2010-05-12', '2010-06-13', '2010-06-29', '2010-07-23', '2010-07-31', '2010-08-08', '2010-08-16', '2010-08-24', '2010-09-09', '2010-09-25', '2010-10-11', '2010-12-06', '2011-04-21', '2011-05-31', '2011-06-16', '2011-06-24', '2011-07-10', '2011-07-18', '2011-07-26', '2011-08-11', '2011-08-19', '2011-08-27', '2011-09-04', '2011-09-28', '2011-10-05', '2011-10-06', '2011-10-22', '2011-10-30', '2011-11-07', '2011-11-15', '2011-12-25', '2012-03-14', '2012-03-30', '2012-04-15', '2012-05-01', '2012-06-02', '2012-06-18', '2012-07-20', '2012-08-05', '2012-09-06', '2012-10-24', '2012-11-09', '2012-12-11', '2013-01-12', '2013-01-28', '2013-03-01', '2013-04-02', '2013-04-05', '2013-04-10', '2013-04-18', '2013-04-26', '2013-05-04', '2013-05-12', '2013-06-05', '2013-06-21', '2013-06-29', '2013-07-07', '2013-07-15', '2013-07-31', '2013-08-08', '2013-08-16', '2013-09-01', '2013-09-09', '2013-09-17', '2013-09-25', '2013-10-03', '2013-10-27', '2013-12-06', '2013-12-30', '2014-02-16', '2014-03-20', '2014-04-05', '2014-05-23', '2014-05-31', '2014-06-08', '2014-06-24', '2014-07-02', '2014-07-18', '2014-07-26', '2014-08-03', '2014-08-11', '2014-08-19', '2014-08-27', '2014-09-04', '2014-09-12', '2014-09-20', '2014-10-14', '2014-10-22', '2014-10-30', '2014-11-23', '2014-12-01', '2014-12-25', '2015-02-03', '2015-02-27', '2015-03-23', '2015-04-08', '2015-04-16', '2015-05-02', '2015-05-26', '2015-06-11', '2015-06-19', '2015-06-27', '2015-07-21', '2015-07-29', '2015-08-06', '2015-08-14', '2015-08-22', '2015-08-30', '2015-09-07', '2015-09-15', '2015-09-23', '2015-10-09', '2015-10-17', '2015-11-10', '2015-11-18', '2015-11-26', '2015-12-04', '2015-12-12', '2015-12-28', '2016-01-13', '2016-01-21', '2016-02-14', '2016-03-17', '2016-03-25', '2016-04-10', '2016-04-18', '2016-04-26', '2016-05-04', '2016-05-28', '2016-06-13', '2016-06-21', '2016-06-29', '2016-07-15', '2016-07-23', '2016-07-31', '2016-08-08', '2016-08-16', '2016-08-24', '2016-09-01', '2016-09-09', '2016-09-17', '2016-10-19', '2016-11-04', '2016-11-12', '2016-11-20', '2016-11-28', '2016-12-14', '2017-04-05', '2017-04-13', '2017-04-29', '2017-05-07', '2017-06-08', '2017-06-16', '2017-06-24', '2017-07-02', '2017-07-10', '2017-07-18', '2017-07-26', '2017-08-03', '2017-08-11', '2017-08-19', '2017-08-27', '2017-09-04', '2017-09-20', '2017-09-28', '2017-10-06', '2017-10-30', '2017-11-15', '2017-12-01', '2017-12-09', '2017-12-17', '2018-01-10', '2018-03-15', '2018-04-24', '2018-05-02', '2018-05-10', '2018-05-18', '2018-05-26', '2018-06-11', '2018-06-19', '2018-06-27', '2018-07-13', '2018-07-29', '2018-08-14', '2018-08-22', '2018-08-30', '2018-09-07', '2018-09-15', '2018-09-23', '2018-10-01', '2018-10-09', '2018-11-02', '2018-11-18', '2018-12-28', '2019-01-05', '2019-02-22', '2019-03-02', '2019-04-11', '2019-04-19', '2019-05-05', '2019-05-21', '2019-05-29', '2019-06-14', '2019-06-22', '2019-06-30', '2019-07-08', '2019-07-16', '2019-08-01', '2019-08-09', '2019-08-17', '2019-08-25', '2019-09-10', '2019-09-17', '2019-09-18', '2019-10-12', '2019-11-13', '2019-12-07', '2020-01-07', '2020-01-08', '2020-01-16', '2020-01-24', '2020-02-17', '2020-03-04', '2020-03-12', '2020-03-20', '2020-04-21', '2020-04-29', '2020-05-30', '2020-05-31']
observationDates = [datetime.datetime.strptime(date, "%Y-%m-%d") for date in observationDates]

monitoringStartDate = datetime.datetime(2018, 1, 1)

mon = bfast.BFASTMonitor(start_monitor=monitoringStartDate, freq=365, k=3, hfrac=0.25, trend=False, level=0.05, backend="python")
mon.fit(indexedRasterTimeSeries, observationDates)

monitoringDates = observationDates[bisect.bisect(observationDates, monitoringStartDate):]
breakDate = monitoringDates[mon.breaks[0, 0]-1]
print(breakDate.timetuple().tm_year,breakDate.timetuple().tm_yday)

2018 122

install.packages("devtools")
devtools::install_github(c("bfast2/strucchange","bfast2/bfast"), quiet = TRUE)

indexedRasterTimeSeries <- c(3141,NA,2153,NA,2791,2795,2109,2139,2077,2200,NA,1260,1301,2496,1429,NA,NA,NA,NA,3122,2446,2428,1647,1266,NA,25,-996,-21,NA,NA,NA,NA,466,1454,NA,NA,2084,2600,NA,2133,1865,NA,596,410,NA,NA,NA,NA,2438,NA,NA,NA,NA,3382,3753,4647,NA,2989,3356,3223,2732,3536,2332,NA,1957,1325,NA,1990,NA,NA,2740,3510,NA,NA,3019,NA,3676,3790,3370,NA,3504,NA,NA,2984,2762,2794,2845,NA,NA,3076,NA,3213,NA,NA,NA,NA,NA,4071,NA,NA,NA,NA,3387,NA,3456,3222,3173,2609,2650,2535,2441,2191,2472,2350,2847,NA,NA,3484,NA,4704,3412,3687,NA,NA,NA,NA,3946,NA,NA,NA,NA,3764,3266,3451,3552,3289,2517,2588,NA,NA,1487,1746,2084,1884,1917,NA,3765,NA,NA,NA,NA,3950,3663,3733,4007,4202,3566,3926,3573,3804,3504,3678,3318,3255,NA,NA,2768,3016,NA,NA,3641,3685,4161,NA,3975,NA,NA,4019,NA,4026,3697,4053,4292,3741,3911,3575,3467,3410,NA,3291,3473,NA,NA,3467,NA,3801,3022,NA,NA,3363,3742,4035,NA,NA,NA,3996,NA,NA,3848,NA,NA,3709,3312,3148,2734,2944,NA,NA,NA,NA,4086,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,NA,4200)
observationDates <- c("2010-02-13","2010-03-09","2010-04-10","2010-04-18","2010-05-12","2010-06-13","2010-06-29","2010-07-23","2010-07-31","2010-08-08","2010-08-16","2010-08-24","2010-09-09","2010-09-25","2010-10-11","2010-12-06","2011-04-21","2011-05-31","2011-06-16","2011-06-24","2011-07-10","2011-07-18","2011-07-26","2011-08-11","2011-08-19","2011-08-27","2011-09-04","2011-09-28","2011-10-05","2011-10-06","2011-10-22","2011-10-30","2011-11-07","2011-11-15","2011-12-25","2012-03-14","2012-03-30","2012-04-15","2012-05-01","2012-06-02","2012-06-18","2012-07-20","2012-08-05","2012-09-06","2012-10-24","2012-11-09","2012-12-11","2013-01-12","2013-01-28","2013-03-01","2013-04-02","2013-04-05","2013-04-10","2013-04-18","2013-04-26","2013-05-04","2013-05-12","2013-06-05","2013-06-21","2013-06-29","2013-07-07","2013-07-15","2013-07-31","2013-08-08","2013-08-16","2013-09-01","2013-09-09","2013-09-17","2013-09-25","2013-10-03","2013-10-27","2013-12-06","2013-12-30","2014-02-16","2014-03-20","2014-04-05","2014-05-23","2014-05-31","2014-06-08","2014-06-24","2014-07-02","2014-07-18","2014-07-26","2014-08-03","2014-08-11","2014-08-19","2014-08-27","2014-09-04","2014-09-12","2014-09-20","2014-10-14","2014-10-22","2014-10-30","2014-11-23","2014-12-01","2014-12-25","2015-02-03","2015-02-27","2015-03-23","2015-04-08","2015-04-16","2015-05-02","2015-05-26","2015-06-11","2015-06-19","2015-06-27","2015-07-21","2015-07-29","2015-08-06","2015-08-14","2015-08-22","2015-08-30","2015-09-07","2015-09-15","2015-09-23","2015-10-09","2015-10-17","2015-11-10","2015-11-18","2015-11-26","2015-12-04","2015-12-12","2015-12-28","2016-01-13","2016-01-21","2016-02-14","2016-03-17","2016-03-25","2016-04-10","2016-04-18","2016-04-26","2016-05-04","2016-05-28","2016-06-13","2016-06-21","2016-06-29","2016-07-15","2016-07-23","2016-07-31","2016-08-08","2016-08-16","2016-08-24","2016-09-01","2016-09-09","2016-09-17","2016-10-19","2016-11-04","2016-11-12","2016-11-20","2016-11-28","2016-12-14","2017-04-05","2017-04-13","2017-04-29","2017-05-07","2017-06-08","2017-06-16","2017-06-24","2017-07-02","2017-07-10","2017-07-18","2017-07-26","2017-08-03","2017-08-11","2017-08-19","2017-08-27","2017-09-04","2017-09-20","2017-09-28","2017-10-06","2017-10-30","2017-11-15","2017-12-01","2017-12-09","2017-12-17","2018-01-10","2018-03-15","2018-04-24","2018-05-02","2018-05-10","2018-05-18","2018-05-26","2018-06-11","2018-06-19","2018-06-27","2018-07-13","2018-07-29","2018-08-14","2018-08-22","2018-08-30","2018-09-07","2018-09-15","2018-09-23","2018-10-01","2018-10-09","2018-11-02","2018-11-18","2018-12-28","2019-01-05","2019-02-22","2019-03-02","2019-04-11","2019-04-19","2019-05-05","2019-05-21","2019-05-29","2019-06-14","2019-06-22","2019-06-30","2019-07-08","2019-07-16","2019-08-01","2019-08-09","2019-08-17","2019-08-25","2019-09-10","2019-09-17","2019-09-18","2019-10-12","2019-11-13","2019-12-07","2020-01-07","2020-01-08","2020-01-16","2020-01-24","2020-02-17","2020-03-04","2020-03-12","2020-03-20","2020-04-21","2020-04-29","2020-05-30","2020-05-31")

monitoringStartDate <- c(2018,1)

pixelTimeSeries <-  bfast::bfastts(indexedRasterTimeSeries, dates = observationDates, type="irregular")
mon <- bfast::bfastmonitor(pixelTimeSeries, monitoringStartDate, formula = response ~ harmon, order = 3, h = 0.25, level = 0.05, history = "all")

year <- floor(mon$breakpoint)
doy <- (mon$breakpoint-year) * 365 + 1
print(c(year,doy))

2018 162

I believe that both runs are run with the same parameters, and they should give the same break date. Am I missing some parameter?

it's it the same vein as #37.

I've corrected some of the bug in my own fork but i'm now seeing another issue. At line

This file should only rely on Python built-in modules like distutils or setuptools.

If not refactored, the lib cannot be installed as in :

pip install -r requirements.txt

I am unable to install bfast (with pip or from GitHub) due to the following issue:

jan@eceo-ws12:~/Documents/bfast$ pip install bfast
Collecting bfast
  Using cached bfast-0.7.tar.gz (306 kB)
  Preparing metadata (setup.py) ... error
  error: subprocess-exited-with-error
  
  × python setup.py egg_info did not run successfully.
  │ exit code: 1
  ╰─> [26 lines of output]
      /tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/monitor/opencl/bfastfinal.py:50290: FutureWarning: In the future `np.bool` will be defined as the corresponding NumPy scalar.
        'numpy_type': np.bool }
      Traceback (most recent call last):
        File "<string>", line 2, in <module>
        File "<pip-setuptools-caller>", line 34, in <module>
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/setup.py", line 16, in <module>
          import bfast
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/__init__.py", line 1, in <module>
          from .models import BFASTMonitor
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/models.py", line 1, in <module>
          from bfast.monitor import BFASTMonitorPython
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/monitor/__init__.py", line 1, in <module>
          from .opencl import BFASTMonitorOpenCL
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/monitor/opencl/__init__.py", line 1, in <module>
          from .base import BFASTMonitorOpenCL
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/monitor/opencl/base.py", line 20, in <module>
          from .bfastfinal import bfastfinal
        File "/tmp/pip-install-z8jn1skv/bfast_c2d55a80d9da45cd8daa95f8a2731d18/bfast/monitor/opencl/bfastfinal.py", line 50290, in <module>
          'numpy_type': np.bool }
                        ^^^^^^^
        File "/home/jan/miniconda3/envs/bfast/lib/python3.11/site-packages/numpy/__init__.py", line 324, in __getattr__
          raise AttributeError(__former_attrs__[attr])
      AttributeError: module 'numpy' has no attribute 'bool'.
      `np.bool` was a deprecated alias for the builtin `bool`. To avoid this error in existing code, use `bool` by itself. Doing this will not modify any behavior and is safe. If you specifically wanted the numpy scalar type, use `np.bool_` here.
      The aliases was originally deprecated in NumPy 1.20; for more details and guidance see the original release note at:
          https://numpy.org/devdocs/release/1.20.0-notes.html#deprecations. Did you mean: 'bool_'?
      [end of output]
  
  note: This error originates from a subprocess, and is likely not a problem with pip.
error: metadata-generation-failed

× Encountered error while generating package metadata.
╰─> See above for output.

note: This is an issue with the package mentioned above, not pip.
hint: See above for details.

By reading your code i realized that you were using the sklearn lib. But according to what I found in the develop branch, you are not using the intersect parameter of the LinearRegression method.

It means that using LInearRegression all-together is useless as the optimisation problem have a closed solution given by

theta = (H^T*H)^-1*H^T*Y

I am proposing a change to the lib using only numpy to compute the coefficients. According to the tests I've performed on my local computer, I get the same results twice as fast (44s with sklearn, 23 without).

It will also reduce the number of dependencies which is alway important.

Hello,

I read this article https://arxiv.org/pdf/1807.01751.pdf but I don't have entirely sure what are the values returned in the means array. Corresponds to the magnitude?

You can use the subset_1 data: https://drive.google.com/drive/folders/1MkvA-MSvZCu_8WHkQ4p3R8Z0Qv_Vq-JC?usp=sharing
And this gist: https://gist.github.com/Carolina710/330b473478dd59188edaaa6a3fd76bb6

My means values are always negative, in Python version is like Max:0, Min:-3.1542912 and OpenCL Max:0, Min:-47.075047

What results could I rely on more accurately? Python or OpenCL?

Thank you all!

I read in the code that the level parameter of BFASTMonitor need to be in [0.95, 0.999]:

__critval_level = np.arange(0.95, 0.999, 0.001)

But I tried to launch it with value .95 I got the following error :

ValueError: ('level can only be one of', array([0.95 , 0.951, 0.952, 0.953, 0.954, 0.955, 0.956, 0.957, 0.958,
       0.959, 0.96 , 0.961, 0.962, 0.963, 0.964, 0.965, 0.966, 0.967,
       0.968, 0.969, 0.97 , 0.971, 0.972, 0.973, 0.974, 0.975, 0.976,
       0.977, 0.978, 0.979, 0.98 , 0.981, 0.982, 0.983, 0.984, 0.985,
       0.986, 0.987, 0.988, 0.989, 0.99 , 0.991, 0.992, 0.993, 0.994,
       0.995, 0.996, 0.997, 0.998, 0.999]))

Which make no sense because the value I've set up is prescisely 0.95.

I looked in your code and found out that the default value is .05... strange. I finally found my answer

def compute_lam(N, hfrac, level, period):
    check(hfrac, period, 1 - level, "max")

    return get_critval(hfrac, period, 1 - level, "max")

You are not testing level but 1-level. Could you please change the error message which is very confusing ?

I am using this implementation of BFAST in the Bfast-gpu module of the SEPAL platform and I wanted to check how they behave on a GPU machine.

I did my tests on a G4 machine (4 GPU, 16 GB of RAM) on a 5 year TS over the Juaboso province (1 370 km² in res 30x30m):

• backend OpenCL: 27 min
• backend Python: 14 min

I was extremely surprised to see the python behaving faster even on GPU machine so I tested it over the test data included in the repository (https://github.com/diku-dk/bfast/blob/master/examples/peru_small_opencl.py):

• backend openCL: 192 s
• backend Python: 43 s

I'm showing you the mean of the 3 best out 10.

My question is simple what is the interest of the OpenCl implementation if it's not faster than the Python one on GPU machines ? Is it normal ?

GPU02 can be used for testing

Hello,

I share "Tile_1" data: https://drive.google.com/drive/folders/15DK5T9hZUaaIBNsfzR3tTp5pFb9XCWt5?usp=sharing

Using this data, an error occurs while processing row 120, in the Python version. This error does not always happen, in this case, it happens when I use start_monitor = datetime(2018, 6, 18), for example. If you try start_monitor = datetime(2016, 1, 1) everything work fine.

The error is the follow:

But running in OpenCL I can obtain results.

You can use this gist: https://gist.github.com/Carolina710/330b473478dd59188edaaa6a3fd76bb6 to work with Tile_1.

Thank you!

Probably nobody needs k>10, but k<3 are needed.

What surprised me even more is that the limitation only exists in the opencl backend.
The only place where k is used in an unsafe way, imo, is on line 88 in the futhark code, but I couldn't understand what is ns, and why it can't be what I assume to be 3, or 5.

Interestingly enough, the corresponding code in the "python" backend, is the same, with the exception that it doesn't differentiate between trend and no trend for the value of k2p2, and there's no limitation for k there.

It would be interesting for me to also understand why the python and opencl implementation calculate sigma slightly differently, perhaps pointing me to some literature would suffice.

But the core issue is still for the opencl implementation to work with k=1 or k=2. If this doesn't make mathematical or technical sense, it should probably be explained in the documentation.