Help me about using ARIMA forecasting rainfall I am currently using the ARIMA provided in R. I use the data as the rainfall time series in QuyNhon (Vietnam) from 2000 to 2017 to forecast rainfall for the next several years.
I wish that the professor could answer for me if I have used the ARIMA model to accurately predict. If the model is wrong, I hope professor can help me rewrite the R code to forecast rainfall more accurately.
I hope that the professor gives me a 95% confidence interval that is correct for rain forecast. And how to use ggplot () to represent this confidence interval.
Link data and code rainfall:
https://drive.google.com/drive/u/1/folders/1MmFeWoUSfXrGNbsL4SiK28FoMAby0gwD
Once again, I thank you very much. I am looking forward to hearing from you soon.
 A: The problem with using a simple arima approach is that model identification is done when it is based upon the tacit assumption that there is no deterministic structure in the data and model parameters are invariant over time.
Your rainfall data (216 monthly values) has two significant seasonal dummies ( October and November) and a few pulses ( 6 of them all in the 4th quarter of the year )  .. one time anomalies and an arima structure (0,0,0)(1,0,0,). You need to integrate both memory (arima) and latent deterministic structure to get a useful/thorough model which has decomposed the observed data to signal and white noise i.e. an error term free of identifiable structure .
Your data  has an ACF  .  A useful model is here  and here  with model summary stats here  with forecasts here for the next 36 months  with an ACF of here suggesting model sufficiency  
The model's residuals are here 
The Actual/Fit and Forecast is here 
You can probably come close to this result by using the TSOUTLIER package in R . I would closely examine the resultant error term to assess it's randomness i.e. the sufficiency of your model.
Hope this helps . Sorry , I do not write R code .
Note that the prediction limits are based upon bootstrapping the model residuals and you can safely truncate the lower values to be => 0.0 .
You asked for 95% prediction limits ..
REVISION OF FORECAST PREDICTION LIMITS:
Using standard computation assuming normality of the errors and no future anomalies this is the Actual/Fit and Forecast  providing more realistic limits ... again simply change the lower limit if it is less than zero to zero.
EDITED IN RESPONSE TO YOUR COMMENT AS TO WHY AUTO.ARIMA developed it's model:
auto.arima can work if there are no latent deterministic structures . Your data has them in a major way. Simply trying a brute-force AIC optimization fails because it uses the ACF while "The correlogram should be calculated from residuals using a model that controls for intervention administration, otherwise the intervention effects are taken to be Gaussian noise, underestimating the actual autoregressive effect." see @Adam0's comment on Interrupted Time Series Analysis - ARIMAX for High Frequency Biological Data?. 
auto.arima incorrectly uses seasonally differences and a seasonal ma term while the correct fix for seasonal non-stationarity (in this case) is to identify and incorporate fixed effects for the two months (Oct and Nov) .  Being unable to identify large pulses and there are a few, it used an ma(1) model memory rather than fixed pulse effects. 
Finally model identification is an iterative , multi-state approach like peeling an onion https://autobox.com/pdfs/ARIMA%20FLOW%20CHART.pdf . 

