<jupyter_start><jupyter_text>eICU Experiments<jupyter_code>import tensorflow as tf
import numpy as np
import h5py
from sklearn.model_selection import train_test_split
import matplotlib.pyplot as plt
import tensorflow_probability as tfp
import sklearn
from sklearn import metrics
import seaborn as sns
import random<jupyter_output><empty_output><jupyter_text>Follow Read-me instruction to download the medical data.After having downloaded the data in '../data/eICU_data.h5', upload the data:<jupyter_code>def get_data(test=True):
    hf = h5py.File('../data/eICU_data.h5', 'r')
    data_total = np.array(hf.get('x'))
    endpoints_total = np.array(hf.get('y'))
    hf.close()
    data_train, data_val, y_train, endpoints_total_val = train_test_split(data_total[:int(len(data_total) * 0.85)],
                                                                          endpoints_total[:int(len(data_total) * 0.85)],
                                                                          test_size=0.20,
                                                                          random_state=42)
    if test:
        data_val = data_total[int(len(data_total) * 0.85):]
        endpoints_total_val = endpoints_total[int(len(data_total) * 0.85):]
    return data_train, data_val, y_train, endpoints_total_val
def batch_generator(data_train, data_val, endpoints_total_val, batch_size, mode="train"):
    while True:
        if mode == "train":
            for i in range(len(data_train) // batch_size):
                time_series = data_train[i * batch_size: (i + 1) * batch_size]
                yield time_series, i
        elif mode == "val":
            for i in range(len(data_val) // batch_size):
                time_series = data_val[i * batch_size: (i + 1) * batch_size]
                time_series_endpoint = endpoints_total_val[i * batch_size: (i + 1) * batch_size]
                yield time_series, time_series_endpoint, i
        else:
            raise ValueError("The mode has to be in {train, val}")<jupyter_output><empty_output><jupyter_text>Insert the name of the job in ex_name:<jupyter_code>ex_name="hyperopt_LSTM_20_16-16_2020-02-17_35a17"<jupyter_output><empty_output><jupyter_text>Get the data:<jupyter_code>batch_size=128
modelpath = "../models/{}/{}".format(ex_name, ex_name)
data_train, data_val, endpoints_total_train, endpoints_total_val = get_data(test=True)<jupyter_output><empty_output><jupyter_text>Create heat-maps, trajectories and probability distributions<jupyter_code>som_dim = [16,16]
latent_dim=10
val_gen = batch_generator(data_train, data_val, endpoints_total_val, 300, mode="val")
num_batches = len(data_val) // 300
tf.reset_default_graph()
with tf.Session() as sess:
    sess.run(tf.global_variables_initializer())
    saver = tf.train.import_meta_graph(modelpath+".meta")
    saver.restore(sess, modelpath)
    graph = tf.get_default_graph()
    k = graph.get_tensor_by_name("k/k:0")
    z_e = graph.get_tensor_by_name("z_e_sample/z_e:0")
    x = graph.get_tensor_by_name("inputs/x:0")
    is_training = graph.get_tensor_by_name("is_training/is_training:0")
    graph = tf.get_default_graph()
    z_e_p = graph.get_tensor_by_name("prediction/next_state/input_lstm:0")
    q = graph.get_tensor_by_name("q/distribution/q:0")
    embeddings = graph.get_tensor_by_name("embeddings/embeddings:0")
    reconstruction = graph.get_tensor_by_name("reconstruction_e/x_hat:0")
    
    print("Evaluation...")
    test_k_all = []
    labels_val_all = []
    z_e_all=[]
    z_q_all = []
    qq = []
    for i in range(num_batches):
            batch_data, batch_labels, ii = next(val_gen)
            f_dic = {x: batch_data}
            test_k_all.extend(sess.run(k, feed_dict=f_dic))
            labels_val_all.extend(batch_labels)
            z_q_all.extend(sess.run(q, feed_dict=f_dic))
            qq.extend(sess.run(q, feed_dict=f_dic))
            z_e_all.extend(sess.run(z_e, feed_dict=f_dic))
    labels_val_all = np.array(labels_val_all)
    k_all = np.array(test_k_all)
    qq = np.array(qq)
    labels_val_all = np.reshape(labels_val_all, (-1, labels_val_all.shape[-1]))
    NMI_24 = metrics.normalized_mutual_info_score(labels_val_all[:, 3], k_all)
    NMI_12 = metrics.normalized_mutual_info_score(labels_val_all[:, 2], k_all)
    NMI_6 = metrics.normalized_mutual_info_score(labels_val_all[:, 1], k_all)
    NMI_1 = metrics.normalized_mutual_info_score(labels_val_all[:, 0], k_all)
    
    embb = sess.run(embeddings, feed_dict={x: data_val})<jupyter_output><empty_output><jupyter_text>Get the labels:<jupyter_code>labels_12 = labels_val_all[:,2]
labels_1 = labels_val_all[:,0]
labels_6 = labels_val_all[:,1]
labels_24 = labels_val_all[:,3]
hosp_disc_1 = labels_val_all[:,4]
hosp_disc_6 = labels_val_all[:,5]
hosp_disc_12 = labels_val_all[:,6]
hosp_disc_24 = labels_val_all[:,7]
u_disc_1 = labels_val_all[:,8]
u_disc_6 = labels_val_all[:,9]
u_disc_12 = labels_val_all[:,10]
u_disc_24 = labels_val_all[:, 11]
labels_1 = labels_1.astype(int)
labels_6 = labels_6.astype(int)
labels_12 = labels_12.astype(int)
labels_24 = labels_24.astype(int)
hosp_disc_12 = hosp_disc_12.astype(int)
hosp_disc_24 = hosp_disc_24.astype(int)
hosp_disc_1 = hosp_disc_1.astype(int)
hosp_disc_6 = hosp_disc_6.astype(int)<jupyter_output><empty_output><jupyter_text>Moran Index<jupyter_code>sd = som_dim[0]*som_dim[1]
mean = np.sum(labels_val_all[:, 0]) / len(labels_val_all[:, 0])
ones = np.ones((len(np.reshape(k_all, (-1)))))
clust_matr1 = np.zeros(som_dim[0]*som_dim[1])
labels= labels_val_all[:, 0]
for i in range(som_dim[0]*som_dim[1]):
    dd = np.sum(ones[np.where(np.reshape(k_all, (-1))==i)])
    if dd == 0:
        s1 = 0
    else:
        s1 = np.sum(labels[np.where(np.reshape(k_all, (-1))==i)]) / np.sum(ones[np.where(np.reshape(k_all, (-1))==i)])
    clust_matr1[i] = s1

k = np.arange(0,sd)
k1 = k // som_dim[0]
k2 = k % som_dim[0]
W = np.zeros((sd,sd))
for i in range(sd):
    for j in range(sd):
        d1 = np.abs((k1[i] - k1[j]))
        d2 = np.abs((k2[i] - k2[j]))
        d1 = min(som_dim[0]-d1, d1)
        d2 = min(som_dim[0]-d2, d2)
        W[i,j] = np.exp(-(d1+d2))
        if i==j:
            W[i,j]=0
M = 0
N_n = 0
for i in range(sd):
    for j in range(sd):
        M += (clust_matr1[i] -mean)*(clust_matr1[j] -mean)* W[i,j]
for i in range(sd):
    N_n += (clust_matr1[i]-mean)**2
W_n = np.sum(W)
I = M * sd / (N_n*W_n)
print(I)<jupyter_output><empty_output><jupyter_text>APACHE score heatmap:<jupyter_code>labels = labels_1
ones = np.ones((len(np.reshape(k_all, (-1)))))
clust_matr1 = np.zeros(som_dim[0]*som_dim[1])
clust_matr2 = np.zeros(som_dim[0]*som_dim[1])
for i in range(som_dim[0]*som_dim[1]):
    s1 = np.sum(labels[np.where(np.reshape(k_all, (-1))==i)]) / np.sum(ones[np.where(np.reshape(k_all, (-1))==i)])
    clust_matr1[i] = s1
clust_matr1 = np.reshape(clust_matr1, (som_dim[0],som_dim[1]))
ax = sns.heatmap(clust_matr1, cmap="YlGnBu")
plt.show()<jupyter_output><empty_output><jupyter_text>Trajectories<jupyter_code>T = []
S = []
for i in range(1000):
    h = np.reshape(u_disc_1, (-1,72))
    if np.max(h[i]) == 1:
        T.append(i)
    else:
        S.append(i)
ind_r = np.random.random_integers(0, 50, 10)
ind_s = np.random.random_integers(0, 50, 10)
T = np.array(T)
S = np.array(S)
a = np.concatenate([S[ind_s], T[ind_r]])
k_all.shape
labels = labels_1
it = 0
k_all = np.reshape(k_all, (-1,72))
fig, ax = plt.subplots(5, 4, figsize=(50,43)) 
ones = np.ones((len(np.reshape(k_all, (-1)))))
clust_matr1 = np.zeros(som_dim[0]*som_dim[1])
clust_matr2 = np.zeros(som_dim[0]*som_dim[1])
for i in range(som_dim[0]*som_dim[1]):
    s1 = np.sum(labels[np.where(np.reshape(k_all, (-1)) == i)]) / np.sum(ones[np.where(np.reshape(k_all, (-1))==i)])
    clust_matr1[i] = s1
clust_matr1 = np.reshape(clust_matr1, (som_dim[0],som_dim[1]))
for t in a:
    #fig, ax = plt.subplots(figsize=(10,7.5)) 
    if it > 9:
        c = "r"
        #print(t)
    else:
        c = "g"
    cc = it % 4
    rr = it // 4
    g = sns.heatmap(clust_matr1, cmap="YlGnBu",ax=ax[rr][cc])
    k_1 = k_all[t] // som_dim[1]
    k_2 = k_all[t] % som_dim[1]
    ax[rr][cc].plot(k_2[:] + 0.5, k_1[:] + 0.5, color=c, linewidth=4)
    ax[rr][cc].scatter(k_2[0] + 0.5, k_1[0] + 0.5, color=c, s=200, label='Start')
    ax[rr][cc].scatter(k_2[1:-1] + 0.5, k_1[1:-1] + 0.5, color=c, linewidth=5, marker='.')
    ax[rr][cc].scatter(k_2[-1] + 0.5, k_1[-1] + 0.5, color=c, s=500, linewidth=4, marker='x', label='End')
    ax[rr][cc].legend(loc=2, prop={'size': 20})
    it +=1
plt.show()<jupyter_output><empty_output><jupyter_text>Probability distribution over trajectory:<jupyter_code>qq.shape
prob_q = np.reshape(qq, (-1, 72, som_dim[0]*som_dim[1])) 
i = np.random.randint(0, 50) #Randomly sampled patient
it = 0
fig, ax = plt.subplots(2, 3, figsize=(50,25))
k_all = np.reshape(k_all, (-1,72))
for t in [0, 17, 40, 57, 64, 71]:
    cc = it % 3
    rr = it // 3
    k_1 = k_all[i] // som_dim[1]
    k_2 = k_all[i] % som_dim[1]
    c = "black"
    g1 = sns.heatmap(np.reshape(prob_q[i, t], (som_dim[0],som_dim[1])), cmap='Reds', alpha=1,  ax=ax[rr][cc])
    ax[rr][cc].plot(k_2[:] + 0.5, k_1[:] + 0.5, color=c, linewidth=6)
    ax[rr][cc].scatter(k_2[0] + 0.5, k_1[0] + 0.5, color=c, s=800, label='Start')
    ax[rr][cc].scatter(k_2[1:-1] + 0.5, k_1[1:-1] + 0.5, color=c, linewidth=10, marker='.')
    ax[rr][cc].scatter(k_2[-1] + 0.5, k_1[-1] + 0.5, color=c, s=1200, linewidth=10, marker='x', label='End')
    ax[rr][cc].legend(loc=2, prop={'size': 30})  
    ax[rr][cc].set_title("Time-step = {}".format(it*14), fontsize=40)
    it +=1
plt.show()<jupyter_output><empty_output><jupyter_text>Unrolling future time-steps and prediction<jupyter_code>def z_dist_flat(z_e, embeddings):
    """Computes the distances between the encodings and the embeddings."""
    emb = np.reshape(embeddings, (som_dim[0]*som_dim[1], -1))
    z = np.reshape(z_e, (z_e.shape[0], 1, latent_dim))
    z = np.tile(z, [1,som_dim[0]*som_dim[1], 1])
    z_dist = np.square(z-emb)
    z_dist_red = np.sum(z_dist, axis=-1)
    return z_dist_red
val_gen = batch_generator(data_train, data_val, endpoints_total_val, 300, mode="val")
tf.reset_default_graph()
num_batches = len(data_val) // 300
latent_dim = 20
num_pred = 6
som = 16*16
max_n_step = 72
with tf.Session() as sess:
    sess.run(tf.global_variables_initializer())
    saver = tf.train.import_meta_graph(modelpath+".meta")
    saver.restore(sess, modelpath)
    graph = tf.get_default_graph()
    k = graph.get_tensor_by_name("k/k:0")
    z_e = graph.get_tensor_by_name("z_e_sample/z_e:0")
    next_z_e = graph.get_tensor_by_name("prediction/next_z_e:0")
    x = graph.get_tensor_by_name("inputs/x:0")
    is_training = graph.get_tensor_by_name("is_training/is_training:0")
    graph = tf.get_default_graph()
    init_1 = graph.get_tensor_by_name("prediction/next_state/init_state:0")
    z_e_p = graph.get_tensor_by_name("prediction/next_state/input_lstm:0")
    state1 = graph.get_tensor_by_name("prediction/next_state/next_state:0")
    q = graph.get_tensor_by_name("q/distribution/q:0")
    embeddings = graph.get_tensor_by_name("embeddings/embeddings:0")
    z_p = graph.get_tensor_by_name('reconstruction_e/decoder/z_e:0')
    reconstruction = graph.get_tensor_by_name("reconstruction_e/x_hat:0")
    
    print("Evaluation...")
    training_dic = {is_training: True, z_e_p: np.zeros((max_n_step * len(data_val), latent_dim)),
                    init_1: np.zeros((2, batch_size, 100)), z_p: np.zeros((max_n_step * len(data_val), latent_dim))}
    k_all = []
    z_e_all=[]
    z_q_all = []
    qq = []
    x_rec = []
    for i in range(num_batches):
            batch_data, batch_labels, ii = next(val_gen)
            f_dic = {x: batch_data}
            k_all.extend(sess.run(k, feed_dict=f_dic))
            z_q_all.extend(sess.run(q, feed_dict=f_dic))
            z_e_all.extend(sess.run(z_e, feed_dict=f_dic))
            qq.extend(sess.run(q, feed_dict=f_dic))
            f_dic.update(training_dic)
            x_rec.extend(sess.run(reconstruction, feed_dict=f_dic))
    z_e_all = np.array(z_e_all)
    k_all = np.array(k_all)
    qq = np.array(qq)
    x_rec = np.array(x_rec)
    z_e_all = z_e_all.reshape((-1, max_n_step, latent_dim))
    k_all = k_all.reshape((-1, max_n_step))
    
    t = 72-num_pred
    
    embeddings = sess.run(embeddings, feed_dict={x: data_val[:, :t, :]})
    embeddings = np.reshape(embeddings,(-1, latent_dim))
    
    z_e_o = z_e_all[:, :t, :]
    k_o = k_all[:, :t]
    k_eval=[]
    next_z_e_o = []
    state1_o =[]
    for i in range(num_batches):
        batch_data, batch_labels, ii = next(val_gen)
        batch_data=batch_data[:, :t, :]
        f_dic = {x: batch_data}
        f_dic.update(training_dic)
        next_z_e_o.extend(sess.run(next_z_e, feed_dict=f_dic))
        if i == 0:
            state1_o = sess.run(state1, feed_dict=f_dic)
        else:
            state1_o = np.concatenate([state1_o, sess.run(state1, feed_dict=f_dic)], axis=1)
    next_z_e_o = np.array(next_z_e_o)
    state1_o = np.array(state1_o)
    
    next_z_e_o_all = np.reshape(next_z_e_o[:, -1, :], (-1,1,latent_dim))
    next_z_e_o = next_z_e_o[:, -1, :]
    k_next = np.argmin(z_dist_flat(next_z_e_o, embeddings), axis=-1)
    k_o = np.concatenate([k_o, np.expand_dims(k_next,1)], axis=1)
    z_e_o = np.concatenate([z_e_o, np.expand_dims(next_z_e_o, 1)], axis=1)
    f_dic = {x: np.zeros((len(data_val),1, 98)), is_training: False, z_e_p: np.zeros((1 * len(data_val), latent_dim)),
             z_p: next_z_e_o, init_1: np.zeros((2, batch_size, 100))}
    x_pred_hat = np.reshape(sess.run(reconstruction, feed_dict=f_dic), (-1, 1, 98))
    
    for i in range(num_pred-1):
        print(i)
        inp = data_val[:1500, (t + i), :]
        f_dic = {x: np.reshape(inp, (inp.shape[0],1,inp.shape[1]))}
        val_dic = {is_training: False, z_e_p: next_z_e_o, init_1: state1_o, z_p: np.zeros((max_n_step * len(inp), latent_dim))}
        f_dic.update(val_dic)
        next_z_e_o = sess.run(next_z_e, feed_dict=f_dic)
        state1_o = sess.run(state1, feed_dict=f_dic)
        next_z_e_o_all = np.concatenate([next_z_e_o_all, next_z_e_o], axis=1)
        k_next = np.argmin(z_dist_flat(next_z_e_o, embeddings), axis=-1)
        k_o = np.concatenate([k_o, np.expand_dims(k_next,1)], axis=1)
        z_e_o = np.concatenate([z_e_o, next_z_e_o], axis=1)
        next_z_e_o = np.reshape(next_z_e_o, (-1, latent_dim))
        f_dic = {x: np.zeros((len(data_val),1, 98)), is_training: False, z_e_p: np.zeros((max_n_step * len(data_val), latent_dim)),
             z_p: next_z_e_o, init_1: np.zeros((2, batch_size, 100))}
        final_x = sess.run(reconstruction, feed_dict=f_dic)
        x_pred_hat = np.concatenate([x_pred_hat, np.reshape(final_x, (-1, 1, 98))], axis = 1)
    
    f_dic = {x: np.zeros((1500,1, 98)), is_training: False, z_e_p: np.zeros((max_n_step * 1500, latent_dim)),
             z_p: z_e_all[:, t-1, :], init_1: np.zeros((2, batch_size, 100))}
    final_x = sess.run(reconstruction, feed_dict=f_dic)
sklearn.metrics.mean_squared_error(np.reshape(x_pred_hat, (-1, 98)), np.reshape(data_val[:1500, -num_pred:], (-1, 98)))<jupyter_output><empty_output><jupyter_text>Accuracy of unrolled state:<jupyter_code>k_true = np.reshape(k_all[:, -num_pred:], (-1))
k_pred = np.reshape(k_o[:, -num_pred:], (-1))
tot = 0
acc = 0
for i in range(len(k_true)):
    tot += 1
    if k_true[i] == k_pred[i]:
        acc += 1
acc = acc / tot
acc<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.1632692681032809, 'fraction_numerical': 0.035074095644858604, 'mean_word_length': 3.7187786259541986, 'pattern_counts': {'":': 2, '<': 49, '<?xml version=': 0, '>': 50, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 12, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '12', 'dir': 'jupyter-structured-clean-dedup', 'id': '8243', 'license': 'MIT', 'n_tokens_mistral': 6880, 'n_tokens_neox': 6604, 'n_words': 1028, 'path': 'notebooks/eICU_experiments.ipynb', 'repo_name': 'jvpoulos/dpsom', 'size': '15453'}