# ---
# jupyter:
#   jupytext:
#     text_representation:
#       extension: .py
#       format_name: light
#       format_version: '1.5'
#       jupytext_version: 1.14.4
#   kernelspec:
#     display_name: BASE
#     language: python
#     name: base
# ---

# # Acquiring Data from open repositories
#
# A crucial step in the work of a computational biologist is not only to analyse data, but acquiring datasets to analyse as well as toy datasets to test out computational methods and algorithms. The internet is full of such open datasets. Sometimes you have to sign up and make a user to get authentication, especially for medical data. This can sometimes be time consuming, so here we will deal with easy access resources, mostly of modest size. Multiple python libraries provide a `dataset` module which makes the effort to fetch online data extremely seamless, with little requirement for preprocessing.
#
# #### Goal of the notebook
# Here you will get familiar with some ways to fetch datasets from online. We do some data exploration on the data just for illustration, but the methods will be covered later.
#
#
# # Useful resources and links
#
# When playing around with algorithms, it can be practical to use relatively small datasets. A good example is the `datasets` submodule of `scikit-learn`. `Nilearn` (library for neuroimaging) also provides a collection of neuroimaging datasets. Many datasets can also be acquired through the competition website [Kaggle](https://www.kaggle.com), in which they describe how to access the data.
#
#
# ### Links
# - [OpenML](https://www.openml.org/search?type=data)
# - [Nilearn datasets](https://nilearn.github.io/modules/reference.html#module-nilearn.datasets)
# - [Sklearn datasets](https://scikit-learn.org/stable/modules/classes.html?highlight=datasets#module-sklearn.datasets)
# - [Kaggle](https://www.kaggle.com/datasets)
# - [MEDNIST]
#
# -  [**Awesomedata**](https://github.com/awesomedata/awesome-public-datasets)
#
#  - We strongly recommend to check out the Awesomedata lists of public datasets, covering topics such as [biology/medicine](https://github.com/awesomedata/awesome-public-datasets#biology) and [neuroscience](https://github.com/awesomedata/awesome-public-datasets#neuroscience)
#
# - [Papers with code](https://paperswithcode.com)
#
# - [SNAP](https://snap.stanford.edu/data/)
#   - Stanford Large Network Dataset Collection  
# - [Open Graph Benchmark (OGB)](https://github.com/snap-stanford/ogb)
#   - Network datasets
# - [Open Neuro](https://openneuro.org/)
# - [Open fMRI](https://openfmri.org/dataset/)

# import basic libraries
import numpy as np
import pandas as pd
from matplotlib import pyplot as plt

# We start with scikit-learn's datasets for testing out ML algorithms. Visit [here](https://scikit-learn.org/stable/modules/classes.html?highlight=datasets#module-sklearn.datasets) for an overview of the datasets.

from sklearn.datasets import fetch_olivetti_faces, fetch_20newsgroups, load_breast_cancer, load_diabetes, load_digits, load_iris

# Load the MNIST dataset (images of hand written digits)

X,y = load_digits(return_X_y=True)

y.shape

X.shape #1797 images, 64 pixels per image


# #### exercise 1.  Make a function `plot` taking an argument (k) to visualize the k'th sample. 
# It is currently flattened, you will need to reshape it. Use `plt.imshow` for plotting. 

# # %load solutions/ex2_1.py
def plot(k):
    plt.imshow(X[k].reshape(8,8), cmap='gray')
    plt.title(f"Number = {y[k]}")
    plt.show()


plot(15); plot(450)

faces = fetch_olivetti_faces()

# #### Exercise 2. Inspect the dataset. How many classes are there? How many samples per class? Also, plot some examples. What do the classes represent? 

# +
# # %load solutions/ex2_2.py

# example solution. 
# You are not expected to make a nice plotting function,
# you can simply call plt.imshow a number of times and observe

print(faces.DESCR) # this shows there are 40 classes, 10 samples per class
print(faces.target) #the targets i.e. classes
print(np.unique(faces.target).shape) # another way to see n_classes

X = faces.images
y = faces.target

fig = plt.figure(figsize=(16,5))
idxs = [0,1,2, 11,12,13, 40,41]
for i,k in enumerate(idxs):
    ax=fig.add_subplot(2,4,i+1)
    ax.imshow(X[k])
    ax.set_title(f"target={y[k]}")
    
# looking at a few plots shows that each target is a single person.
# -

# Once you have made yourself familiar with the dataset you can do some data exploration with unsupervised methods, like below. The next few lines of code are simply for illustration, don't worry about the code (we will cover unsupervised methods in submodule F).

from sklearn.decomposition import randomized_svd

X = faces.data

n_dim = 3
u, s, v = randomized_svd(X, n_dim)


# Now we have factorized the images into their constituent parts. The code below displays the various components isolated one by one.

def show_ims(ims):
    fig = plt.figure(figsize=(16,10))
    idxs = [0,1,2, 11,12,13, 40,41,42, 101,101,103]
    for i,k in enumerate(idxs):
        ax=fig.add_subplot(3,4,i+1)
        ax.imshow(ims[k])
        ax.set_title(f"target={y[k]}")


for i in range(n_dim):
    my_s = np.zeros(s.shape[0])
    my_s[i] = s[i]
    recon = u@np.diag(my_s)@v
    recon = recon.reshape(400,64,64)
    show_ims(recon)

# Are you able to see what the components represent? It at least looks like the second component signifies the lightning  (the light direction), the third highlights eyebrows and facial chin shape.

from sklearn.manifold import TSNE

tsne = TSNE(init='pca', random_state=0)
trans = tsne.fit_transform(X)

# +
m = 8*10 # choose 4 people

plt.figure(figsize=(16,10))
xs, ys = trans[:m,0], trans[:m,1]
plt.scatter(xs, ys, c=y[:m], cmap='rainbow')

for i,v in enumerate(zip(xs,ys, y[:m])):
    xx,yy,s = v 
    #plt.text(xx,yy,s) #class
    plt.text(xx,yy,i) #index
# -

# Many people seem to have multiple subclusters. What is the difference between those clusters? (e.g. 68,62,65 versus the other 60's)

# +
ims = faces.images

idxs = [68,62,65,66,60,64,63]
#idxs = [9,4,1, 5,3]
for k in idxs:
    plt.imshow(ims[k], cmap='gray')
    plt.show()


# -

def show(im):
    return plt.imshow(im, cmap='gray')


import pandas as pd
df= pd.read_csv('data/archive/covid_impact_on_airport_traffic.csv')

df.shape

df.describe()

df.head()

df.Country.unique()

df.ISO_3166_2.unique()

df.AggregationMethod.unique()

# Here we will look at [OpenML](https://www.openml.org/) - a repository of open datasets free to explore data and test methods.
#
# ### Fetching an OpenML dataset
#
# We need to pass in an ID to access, as follows:

from sklearn.datasets import fetch_openml

# OpenML contains all sorts of datatypes. By browsing the website we found a electroencephalography (EEG) dataset to explore: 

data_id = 1471 #this was found by browsing OpenML
dataset = fetch_openml(data_id=data_id, as_frame=True)

dir(dataset)

dataset.url

type(dataset)

print(dataset.DESCR)

original_names = ['AF3',
 'F7',
 'F3',
 'FC5',
 'T7',
 'P',
 'O1',
 'O2',
 'P8',
 'T8',
 'FC6',
 'F4',
 'F8',
 'AF4']

dataset.feature_names

df = dataset.frame

df.head()

df.shape[0] / 117
# 128 frames per second

df = dataset.frame
y = df.Class
#df.drop(columns='Class', inplace=True)

df.dtypes

# +
#def summary(s):
#    print(s.max(), s.min(), s.mean(), s.std())
#    print()
#    
#for col in df.columns[:-1]:
#    column = df.loc[:,col]
#    summary(column)
# -

df.plot()

# From the plot we can quickly identify a bunch of huge outliers, making the plot look completely uselss. We assume these are artifacts, and remove them.

df2 = df.iloc[:,:-1].clip_upper(6000)
df2.plot()

# Now we see better what is going on. Lets just remove the frames corresponding to those outliers

frames = np.nonzero(np.any(df.iloc[:,:-1].values>5000, axis=1))[0]
frames

df.drop(index=frames, inplace=True)

df.plot(figsize=(16,8))
plt.legend(labels=original_names)



df.columns

# ### Do some modelling of the data

from sklearn.linear_model import LogisticRegression

lasso = LogisticRegression(penalty='l2')

X = df.values[:,:-1]
y = df.Class
y = y.astype(np.int) - 1 # map to 0,1

print(X.shape)
print(y.shape)

lasso.fit(X,y)

comp = (lasso.predict(X) == y).values

np.sum(comp.astype(np.int))/y.shape[0] # shitty accuracy

lasso.coef_[0].shape

names = dataset.feature_names

original_names

# +
coef = lasso.coef_[0]
plt.barh(range(coef.shape[0]), coef)
plt.yticks(ticks=range(14),labels=original_names)

plt.show()
# -

# Interpreting the coeficients: we naturally tend to read the magnitude of the coefficients as feature importance. That is a fair interpretation, but currently we did not scale our features to a comparable range prior to fittting the model, so we cannot draw that conclusion.

# ### Extra exercise. Go to [OpenML](https://openml.org) and use the search function (or just look around) to find any dataset that interest you. Load it using the above methodology, and try to do anything you can to understand the datatype, visualize it etc.

# +
### YOUR CODE HERE
# -



StarcoderdataJupyterScriptsDedupFiltered