import matplotlib.pyplot as plt
import numpy as np
import pints
import pints.toy

"""
Explanation: Fitting a logistic model to German credit data
This notebook explains how to run the toy logistic regression model example using the German credit data from [1]. In this example, we have predictors for 1000 individuals and an outcome variable indicating whether or not each individual should be given credit.
[1] "UCI machine learning repository", 2010. A. Frank and A. Asuncion. https://archive.ics.uci.edu/ml/datasets/statlog+(german+credit+data)
End of explanation
"""


logpdf = pints.toy.GermanCreditLogPDF(download=True)

"""
Explanation: To run this example, we need to first get the data from [1] and process it so we have dichtonomous $y\in{-1,1}$ outputs and the matrix of predictors has been standardised. In addition, we also add a column of 1s corresponding to a constant term in the regression.
If you are connected to the internet, by instantiating with x=None, Pints will fetch the data from the repo for you. If, instead, you have local copies of the x and y matrices, these can be supplied as arguments.
End of explanation
"""


def pca(X):
  # Data matrix X, assumes 0-centered
  n, m = X.shape
  # Compute covariance matrix
  C = np.dot(X.T, X) / (n-1)
  # Eigen decomposition
  eigen_vals, eigen_vecs = np.linalg.eig(C)
  # Project X onto PC space
  X_pca = np.dot(X, eigen_vecs)
  return X_pca

x, y = logpdf.data()
scores = pca(x)

# colour individual points by whether or not to recommend them credit
plt.scatter(scores[:, 0], scores[:, 1], c=y)
plt.xlabel('PCA 1')
plt.ylabel('PCA 2')
plt.show()

"""
Explanation: Let's look at the data: x is a matrix of predictors and y is a vector of credit recommendations for 1000 individuals. 
Specifically, let's look at the PCA scores and plot the first two dimensions against one another. Here, we see that the two groups overlap substantially, but that there neverless some separation along the first PCA component.
End of explanation
"""


xs = [
    np.random.uniform(0, 1, size=(logpdf.n_parameters())),
    np.random.uniform(0, 1, size=(logpdf.n_parameters())),
    np.random.uniform(0, 1, size=(logpdf.n_parameters())),
]

mcmc = pints.MCMCController(logpdf, len(xs), xs, method=pints.HamiltonianMCMC)
mcmc.set_max_iterations(200)

# Set up modest logging
mcmc.set_log_to_screen(True)
mcmc.set_log_interval(10)

for sampler in mcmc.samplers():
    sampler.set_leapfrog_step_size(0.2)
    sampler.set_leapfrog_steps(10)

# Run!
print('Running...')
chains = mcmc.run()
print('Done!')

"""
Explanation: Now we run HMC to fit the parameters of the model.
End of explanation
"""


results = pints.MCMCSummary(chains=chains, time=mcmc.time())
print(results)

"""
Explanation: HMC is quite efficient here at sampling from the posterior distribution.
End of explanation
"""


