IDEMSInternational / R-Instat

A statistics software package powered by R
http://r-instat.org/
GNU General Public License v3.0
38 stars 102 forks source link

More on the SPI dialogue #8244

Open WezzieBanda opened 1 year ago

WezzieBanda commented 1 year ago

I am getting a bug in the SPI dialogue. Initially, I was using Malawi data which has 4 stations, I then used Nigeria-Samaru 56 data from the library and I have noted that I am getting the same bug.

I went to: Cimatic>Prepare>SPI

IMG_20230328_140957 IMG_20230328_141014

rdstern commented 1 year ago

@Vitalis95 or an alternative person @N-thony I second the problem with this dialog. However this seems to be caused by the new code for the output window, so it should be easy to fix?

I tried with the single station, wichita, provided in the SPEI package.

image

In the old Version, 0.7.6 it works and here is the R-code:

# Dialog: SPI/SPEI

wichita <- data_book$get_data_frame(data_name="wichita", use_current_filter=FALSE)
data_ts <- spei_input(year="YEAR", month="MONTH", element="PRCP", data=wichita)
spi_mod1 <- SPEI::spi(data=data_ts, scale=1, kernel=list(type='gaussian', shift=0), na.rm=TRUE)
data_book$add_model(model_name="spi_mod1", model=spi_mod1, data_name="wichita")
spi <- spei_output(x=data_book$get_models(data_name="wichita", model_name="spi_mod1"), year="YEAR", month="MONTH", data=wichita)
data_book$add_columns_to_data(data_name="wichita", col_name="spi", col_data=spi, before=FALSE)

rm(list=c("spi", "spi_mod1", "data_ts", "wichita")) 

And this works fine,

In the new version it gives the error above and the R code is now:

# Dialog: SPI/SPEI

wichita <- data_book$get_data_frame(data_name="wichita", use_current_filter=FALSE)
data_ts <- spei_input(year="YEAR", month="MONTH", element="PRCP", data=wichita)
spi_mod3 <- SPEI::spi(data=data_ts, scale=1, kernel=list(type='gaussian', shift=0), na.rm=TRUE)
data_book$add_object(data_name="wichita", object_name="spi_mod3", object_type_label="model", object_format="text", object=spi_mod3)
wichita <- data_book$get_data_frame(data_name="wichita", use_current_filter=FALSE)
spi <- spei_output(x=data_book$get_object_data(data_name="wichita", object_name="spi_mod3", as_file=TRUE), year="YEAR", month="MONTH", data=wichita)
data_book$add_columns_to_data(data_name="wichita", col_name="spi", col_data=spi, before=FALSE)

rm(list=c("spi", "spi_mod3", "data_ts", "wichita"))

I assume that should be easy to fix. However, the method in the package now has a plot method and so can also produce a ggplot of the resulting index.

image

So could there also be a checkbox on the dialog to give, and perhaps save the plot.

I just added plot(spi_mod1) to the Version 0.7.6 code to give that plot. However, the dialog should work with multiple data sets (more than one station), so the ggplot code needs to be a bit more general than that.

However I now think we sghould use our own ggplot dialog a bit like Danny has done for the PICSA graphs. Here are the ideas.

So perhaps leave this, or for now just do the default plot. Then do our own later. In the manual it says:

# Modding the plot
# Since plot.spei() returns a ggplot object, it is possible to add or tweak
# parts of the plot.
require(ggplot2)
plot(spei(wichita[, "BAL"], 12)) +
ggtitle("SPEI1 at Wichita") +
scale_fill_manual(values = c("blue", "red")) + # classic SPEI look
scale_color_manual(values = c("blue", "red")) + # classic SPEI look
theme_classic() +
theme(legend.position = "bottom")

So I assume you could easily include a facet option too?

rdstern commented 1 year ago

@Vitalis95 the bug is now fixed and I hope that version can be merged for the forthcoming release.

However, the graph option is not great. Let's not do anything to that now. I would like to test the current implementation with Chiara on the workshop in CIMH and then I hope she and they will agree with the following implementation.

With the current version the graph has the same colour above and below the axis. Instead I would like to do the following work.

a) Calculate a new variable called spi_plus with the calculator, namely spi > 0.
b) Then use the bar chart dialog seconf button, and by spi_plus.

That gives a nice plot. And we can have facets if there are multiple stations.

Then discuss the proposed changes below, namely: c) Delete the plot option from the dialog.
d) Add a checkbox to Save Plus (I am not sure what to call this? If checked, then it saves the calculated value above.

Then we could add plots to the dialog, but it is easy enough to use the bar chart and I also wonder about having plots in the Climatic > Describe menu later. In that case, all this does is sort of "what it says on the tin" namely produces the summary.

I notice the spi object has plot implemented. So we could also see how to implement that as a standard way of dealing with objects.

rdstern commented 10 months ago

@fran2or and @jkmusyoka I suggest we could do more on the SPI dialog and could perhaps start by formulating 1 or 2 mini-projects for the AIMS students?

a) The SPEI package has been substantially updated in 2023. I suggest we should provide an improved dialog for R-Instat. I also suggest we move it. In the climatic menu it is currently in the Climatic > Prepare menu. I suggest we move it to the Climatic Describe menu. b) The github site for SPEI includes a specifice problem on the use of SPEI with daily data. It is here.

I replied to her and she has supplied the data she used. It is interesting data, because it is for 316 stations.

Can the group use her data and solve her problem, or at least investigate her problem?

I searched and found she used the dialog in her MSc thesis. This looks to be an impressive thesis and is here.

image