How to add confidence to model's prediction?
It is common to say the error term follows a standard Guassian distribution. If you assume that to be true, then your squared errors follow a Chi-squared distribution:


  In probability theory and statistics, the chi-squared distribution (also chi-square or χ2-distribution) with k degrees of freedom is the distribution of a sum of the squares of k independent standard normal random variables.


Have a look here for some ideas about how to implement a quasi-confidence metric, based on your (mean squared) errors. It assumes the errors follow a chi-aquared distribution and then uses the normalised RMSE to define a set of confidence boundaries for a given confidence level,  $\alpha$, as follows:

$$
\left[\sqrt{\frac{n}{\chi_{1-\frac{\alpha}{2},n}^{2}}}RMSE,\sqrt{\frac{n}{\chi_{\frac{\alpha}{2},n}^{2}}}RMSE\right]
$$

See the link for the steps involved. Here is the coded simulation taken from that post, with some added comments (requires python 3):

from scipy import stats
import numpy as np

s = 3                           # a constant to scale the random distribution
n = 4                           # number of samples/states per prediction
alpha = 0.05                    # confidence interval

# distribution with confidence intervals ɑ = 0.05
c1, c2 = stats.chi2.ppf([alpha/2, 1-alpha/2], n)

# we will take this many samples (this pre-allocates the y-vector)
y = np.zeros(50000)

# Loop over each sample and record the result mean sample
# This would be your prediction vector - here it is random noise
for i in range(len(y)):
    y[i] = np.sqrt(np.mean((np.random.randn(n)*s)**2))

# Use the chi-squared distributed confidence intervals to see when predictions fall
# finds percentage of samples that are inside the confidence interval
conf = mean((sqrt(n/c2)*y &lt; s) &amp; (sqrt(n/c1)*y &gt; s))

print("1-alpha={:2f}".format(conf))


Here is another answer on CrossValidated, which gives more information around the area. 



Additionally, if you assume your predictions lie within a Gaussian distribution, you could use the variance of your predictions as the confidence (welcome to Bayesian learning!). 

There are packages that will help you do this, such as BayesOptimization. There are lots of examples on that webpage. Essentially, you will be able to make predictions and automatically get robust estimates for condifence... and some cool plots to show where your model is quite sure, and where it isn't: