<jupyter_start><jupyter_text>In this project, we will predict RBP-binding RNA sequences using a neural network model based on CNN (<NAME>, et al. "Predicting the sequence specificities of DNA- and RNA-binding proteins by deep learning", Nature Biotechnology, 33, pages831–838 (2015)). 1) Build CNN models for predicting the RBP-binding RNA sequences (ALKBH5 and ELAVL1).2) Train your model using the training dataset (and optimize your model using validation set).3) Report the ROC-AUC score of the test dataset using the trained model.4) Use google colab (or whatever you want)5) Upload  jupyter notebook (".ipynb" file) or link to your google colab code Due data: Oct/11th 23:59pmGrading criteria (10pt total):- It runs anyway... --> 5pt- Properly optimize the models --> 2pt- Good performance --> 2pt- Nice and proper comments --> 1pt- It's a beautiful code! --> +3pt Project 2에서는 RBP-binding RNA sequence를 예측하는 것인데, 이를 확장해서 binding motif도 역시 구하는 extended project를 수행함프로젝트 수행 수준에 따라, 자네의 최종 점수는 중간고사의 class 평균점수 + alpha (물론 여기서 alpha는 negative number도 가능하다)  Prediction of ELAVL1 binding sequencesIn this notebook, we **predict ELAVL1 binding RNA sequences with a CNN model.**Then as an extention, we also **analyze learned CNN filters to find local binding motif sequences.**To effectively detect binding motifs of different lengths, we used an adapted version of **inception module** used for GoogleNet in image classification task.The notebook is written as follows:**Basic project**- Module / data import- Train / validation set split- CNN architecture design- Hyperparameter setting & model optimization- Model evaluation on test dataset with ROC-AUC score**Extended project**- **Analysis of learned filters for local RNA sequence motifs**References:[1] <NAME>, et al. "Predicting the sequence specificities of DNA- and RNA-binding proteins by deep learning", Nature Biotechnology, 33, pages831–838 (2015)[2] Szegedy et al., "Going Deeper with Convolutions", CVPR (2015)[3] Chung et al., "Prediction of binding property of RNA-binding proteins using multi-sized filters and multi-model deep convolutional network", PLOS ONE[4] Zhang et al., "A deep learning framework for modeling structuralfeatures of RNA-binding protein targets", Nucleic Acids Research Module / data import<jupyter_code>#### project2_pseudo_code.py

# import
import os
import numpy as np
import random

import keras
import tensorflow as tf
###################
import tensorflow.keras as keras
from keras import Model
from keras.backend import one_hot
from keras.layers import Activation, Dense, Dropout, Input, BatchNormalization
from keras.layers import Conv1D, Concatenate, MaxPooling1D, AveragePooling1D, Flatten
from keras.regularizers import l2
from keras.optimizers import adam
from keras.utils import np_utils
###################
# import my google drive
from google.colab import drive
drive.mount('/content/gdrive')
# read fasta files  
fasta_file_path='gdrive/My Drive/data sample/'
rbp_name = 'ELAVL1'
print('list of files', os.listdir(fasta_file_path))
# loading training set
train_fasta = list()
train_label = list()
for single_file in [x for x in os.listdir(fasta_file_path) if rbp_name in x and 'train' in x]:
  print('Processing file...', single_file)
  with open(fasta_file_path + single_file) as f:
    for line in f.readlines():
      # get fasta sequence
      if '>' in line:
        continue
      else:
        train_fasta.append(line.strip())
      # get positive negative label
      if 'positives' in single_file:
        train_label.append(1)
      else:
        train_label.append(0)
# convert sequence file to one-hot encoding representation
# function: fasta to onehot representation
def convert2onehot(sequence_list):
  map = {
      'A':[1,0,0,0],
      'U':[0,1,0,0],
      'T':[0,1,0,0],
      'G':[0,0,1,0],
      'C':[0,0,0,1]
  }
  
  onehot = []
  for single_sequence in sequence_list:
    single_onehot = []
    for x in single_sequence:
      single_onehot.append(map[x.upper()])
    onehot.append(single_onehot)
    
  return np.asarray(onehot, dtype=np.float32)

data_input = convert2onehot(train_fasta)
data_label = keras.utils.to_categorical(train_label, 2)<jupyter_output><empty_output><jupyter_text>Train / validation set split80% training set / 20% validation set<jupyter_code># Install sklearn package
!pip install imblearn
# split training set into training set and validation set
# random shuffling of training data
###################
# random shuffle indices for training data
arr = np.arange(data_input.shape[0])
np.random.shuffle(arr)

# apply same shuffling to data and label
data_input = data_input[arr,:,:]
data_label = data_label[arr,:]

# split training and validation dataset from raw dataset (9:1)
from sklearn.model_selection import train_test_split
train_input, validation_input, train_label, validation_label = train_test_split(data_input, data_label, test_size=0.1, random_state=1)
###################
print('Dataset preparation done... train_input, train_label, validation_input, validation_label')
print('Size of each set...', train_input.shape, train_label.shape, validation_input.shape, validation_label.shape)<jupyter_output><empty_output><jupyter_text>CNN architecture design**Design purpose: Build a model that predicts from & learns various lengths of sequence motifs.**Inception module from GoogleNet does this (concatenation of different sized filters, [2]).![inception_module](https://d2mxuefqeaa7sj.cloudfront.net/s_8C760A111A4204FB24FFC30E04E069BD755C4EEFD62ACBA4B54BBA2A78E13E8C_1490879611424_inception_module.png)We can adapt this architecture to sequence model simply by changing 2d filters to 1d.At the same time, we want to analyze learned filters in 1st inception module to look for RBP binding sequence motifs.From the requirements, we designed the model as follows:- Use simplified inception modules without 1-convolution sublayers- Instead of small (1-, 3-, 5-) convolutional filters, use large ones (8-, 16-, 32-) for later analysis of filters in the first inception module [3]- 3 inception modules -> Input size reduction by max-pooling -> FC layers for classificationOptimization / regularization techniquesBelow are some additional techniques we applied for optimization and regularization:- Adam optimizer- L2 weight regularization- Dropout before prediction layer- Early stoppingModel training curve![training_curve](https://drive.google.com/uc?export=download&id=1xgGiz_9_ECWh5_p8EaiL_V-wRw2jEeJN)<jupyter_code>train_mode = True
# model building
###################
# hyperparameters
lr = 3e-4
reg = 1e-2
p = 1/2

# design purpose: build a model that predicts from / learns various lengths of sequence motifs.
# GoogleNet inception module-based CNN architecture
# use modified inception module instead of plain convolutional layers

# inception module: detects motifs with different kernel sizes
# concatenation of activation maps from filters of variying lengths

# use 1d convolution layers, starting from 4-channel input initially
inputs = Input(shape=(200,4,))

# inception module 1
x1 = Conv1D(2, 8, strides=1, padding='same', data_format = 'channels_last')(inputs)
x2 = Conv1D(2, 16, strides=1, padding='same', data_format = 'channels_last')(inputs)
x3 = Conv1D(2, 32, strides=1, padding='same', data_format = 'channels_last')(inputs)
# concatenate motif scans along channel dimension
x = Concatenate(axis=2)([x1, x2, x3])
# ReLU activation
x = Activation('relu')(x)
# max-pooling
x = MaxPooling1D(pool_size=2, padding='valid', data_format='channels_last')(x)

# inception module 2
x1 = Conv1D(2, 8, strides=1, padding='same', data_format = 'channels_last')(x)
x2 = Conv1D(2, 16, strides=1, padding='same', data_format = 'channels_last')(x)
x3 = Conv1D(2, 32, strides=1, padding='same', data_format = 'channels_last')(x)
# concatenate motif scans along channel dimension
x = Concatenate(axis=2)([x1, x2, x3])
# ReLU activation
x = Activation('relu')(x)
# max-pooling
x = MaxPooling1D(pool_size=2, padding='valid', data_format='channels_last')(x)

# inception module 3
x1 = Conv1D(2, 8, strides=1, padding='same', data_format = 'channels_last')(x)
x2 = Conv1D(2, 16, strides=1, padding='same', data_format = 'channels_last')(x)
x3 = Conv1D(2, 32, strides=1, padding='same', data_format = 'channels_last')(x)
# concatenate motif scans along channel dimension
x = Concatenate(axis=2)([x1, x2, x3])
# ReLU activation
x = Activation('relu')(x)
# max-pooling
x = MaxPooling1D(pool_size=2, padding='valid', data_format='channels_last')(x)

# reducing input size with 1-convolution and max pooling
x = Conv1D(6, 1, strides=1, data_format='channels_last')(x)
x = MaxPooling1D(pool_size=5, padding='valid', data_format='channels_last')(x)
x = Activation('relu')(x)
x = Conv1D(6, 1, strides=1, data_format='channels_last')(x)
x = MaxPooling1D(pool_size=5, padding='valid', data_format='channels_last')(x)
x = Activation('relu')(x)

# flattening before FC layers
x = Flatten()(x)
# plain FC layers
x = Dense(64, kernel_regularizer=l2(reg))(x)
x = Activation('relu')(x)
x = Dropout(rate=p)(x)
# prediction layer
predictions = Dense(2, activation='softmax')(x)

# model compilation
model = Model(inputs=inputs, outputs=predictions)
model.compile(optimizer=adam(lr=lr,decay=1e-5),
                 loss='categorical_crossentropy',
                 metrics=['accuracy'])
model.summary() # print summary of model
###################

# model training
###################
if train_mode:
  # if training mode, train
  # early stopping with validation dataset
  from keras.callbacks import EarlyStopping
  from keras.callbacks import ModelCheckpoint

  es = EarlyStopping(monitor='val_loss', mode='min', verbose=1, patience=5)
  mc = ModelCheckpoint(fasta_file_path + 'ELAVL1_best.h5', monitor='val_acc', mode='max', verbose=1, save_best_only=True)

  # training the model using training and validation set
  epochs = 50
  # train AR-classification task
  hist = model.fit(train_input, train_label,
                   validation_data=(validation_input, validation_label),
                   batch_size=8, epochs=epochs, verbose=0, shuffle=True, callbacks=[es, mc])
  
  # visualize training curve
  import matplotlib.pyplot as plt
  fig, axs = plt.subplots(1,2, constrained_layout=True)

  axs[0].plot(hist.history['loss'], label='train loss')
  axs[0].plot(hist.history['val_loss'], label='val loss')
  axs[0].legend(loc='upper left')
  axs[0].set_title('loss')

  axs[1].plot(hist.history['acc'], label='train acc')
  axs[1].plot(hist.history['val_acc'], label='val acc')
  axs[1].legend(loc='upper left')
  axs[1].set_title('acc')

else:
  # if not training mode, load pre-trained weights
  model.load_weights(fasta_file_path + 'ELAVL1_best.h5')
  print("Model weights loaded")
###################<jupyter_output><empty_output><jupyter_text>Model evaluation on test dataset with ROC-AUC scoreWith trained model, we report the performance on test set as follows:Classification accuracy: 93%AUC-ROC score: 0.97![roc_curve](https://drive.google.com/uc?export=download&id=1igHW3s2RsPE3m_y0ZA2ILKZ_RxJj1jON)<jupyter_code># evaluation: calculate AUC score for the test set
# loading test sets 
###################
# loading training set
test_fasta = list()
test_label = list()
for single_file in [x for x in os.listdir(fasta_file_path) if rbp_name in x and 'test' in x]:
  print('Processing file...', single_file)
  with open(fasta_file_path + single_file) as f:
    for line in f.readlines():
      # get fasta sequence
      if '>' in line:
        continue
      else:
        test_fasta.append(line.strip())
      # get positive negative label
      if 'positives' in single_file:
        test_label.append(1)
      else:
        test_label.append(0)
###################

test_input = convert2onehot(test_fasta)
test_label = keras.utils.to_categorical(test_label, 2)
print('TEST set prepared...', test_input.shape, test_label.shape)

# calculate AUC score for test sets
###################
# evaluate the prediction model using test dataset.
from sklearn.metrics import classification_report, roc_curve, auc
import matplotlib.pyplot as plt

test_loss, test_acc = model.evaluate(test_input, test_label)
test_pred = model.predict(test_input)

# 0: negative, 1: positive
target_names = ['0','1']
print("classification report")
print(classification_report(test_label[:,1]==1, test_pred[:,1]>0.5, target_names=target_names))

# compute and draw ROC curve
fpr, tpr, thresholds = roc_curve(test_label[:,1], test_pred[:,1])
roc_auc = auc(fpr, tpr)

fig = plt.figure()
plt.plot(fpr, tpr, "b")
plt.plot([0, 0], [0, 1], "k--")
plt.plot([0, 1], [1, 1], "k--")
plt.plot([0, 1], [0, 1], "k--")
plt.title("AUC = %0.2f" % roc_auc)
plt.xlabel("FP rate (1 - specificity)")
plt.ylabel("TP rate (sensitivity)")
###################<jupyter_output><empty_output><jupyter_text>Extended project: Analysis of learned filters for local RNA sequence motifsIn this section, we  (1) analyze learned filters of various lengths in our inception-like CNN,  (2) try to find RNA sequence motifs important for RBP binding and  (3) assess their biological plausibility.Convolution filters, or motif detectors learned from RBP data, are shown below as 4 x length colormaps (each columns were mean-subtracted before visualization).We can see that, for example, the first 8-length filter detects U-rich local region and the first 16-length filter specifically detects "CUCCUCCAUGG...".![f8](https://drive.google.com/uc?export=download&id=1fnsSsxMJ0UATQquwDsx1zaQlN0RY3UST)![f16](https://drive.google.com/uc?export=download&id=1TXNbLotXTMkSEFLTOYWjTQ9ezySponHH)![f32](https://drive.google.com/uc?export=download&id=1uKDrw86-gsX_IrnYZOVb_ZFPcJZiV2lz)<jupyter_code># load convolution filters from learned model
conv1 = model.layers[1].get_weights()
conv2 = model.layers[2].get_weights()
conv3 = model.layers[3].get_weights()
# load w and b
w1 = conv1[0]; b1 = conv1[1]
w2 = conv2[0]; b2 = conv2[1]
w3 = conv3[0]; b3 = conv3[1]

# function for visualization of learned filters
def vis_filter(length, w):
  fig = plt.figure()
  plt.title("filter length = "+str(length))
  plt.axis('off')
  for ii in range(2):
    ax = fig.add_subplot(2,1,ii+1)
    temp_img = w[:,:,ii].transpose() # filter image
    ax.imshow(temp_img-np.mean(temp_img,0),cmap='bwr') # mean subtraction for each position
    ax.set_xticks([],[])
    plt.yticks((0,1,2,3), ('A','U','C','G'))

# plot w of 8-length filter as a heatmap
vis_filter(8, w1)
# plot w of 16-length filter as a heatmap
vis_filter(16, w2)
# plot w of 32-length filter as a heatmap
vis_filter(32, w3)<jupyter_output><empty_output><jupyter_text>(source code for sequence logo)modified from: https://github.com/azofeifa/PWM_logo/tree/master/src<jupyter_code>from math import log
import numpy as np
import math
import matplotlib.pyplot as plt

def despine(ax,right=True, top=True, left=False,bottom=False):
    if left:
        ax.spines['left'].set_visible(False)
    if bottom:
        ax.spines['bottom'].set_visible(False)
    if right:
        ax.spines['right'].set_visible(False)
    if top:
        ax.spines['top'].set_visible(False)

def draw_A(width, height, x, y, ax, color="blue",lw=2.0,alpha=0.5):
    ax.plot([x, x+width*0.5], [y, y+height],lw=lw, color=color ,alpha=alpha)
    ax.plot([x+width*0.5, x+width], [y+height, y],lw=lw, color=color ,alpha=alpha)
    ax.plot([x+width*0.25,x+width*0.75], [height*0.5, height*0.5], color=color, lw=lw,alpha=alpha)
    return y+height

def draw_U(width, height,x,y, ax, color="red",lw=2.0,alpha=0.5):
    ax.plot([x, x], [y+height/2, y+height],lw=lw, color=color,alpha=alpha)
    ax.plot([x+width, x+width], [y+height/2, y+height],lw=lw, color=color,alpha=alpha)
    x += width/2
    y += height/2
    ts = np.linspace(math.pi, math.pi*2, 100)
    xs = [math.cos(t)*width/2+x for t in ts]
    ys = [math.sin(t)*height/2+y for t in ts]
    ax.plot(xs,ys,lw=lw, color=color,alpha=alpha)
    return y + height/2

def draw_C(width, height,x,y, ax, color="green",lw=2.0,alpha=0.5):
    x += width/2
    y += height/2
    ts = np.linspace(3.85*math.pi/2., math.pi/6., 100)
    xs = [math.cos(t)*width/2+x for t in ts]
    ys = [math.sin(t)*height/2+y for t in ts]
    ax.plot(xs,ys,lw=lw, color=color,alpha=alpha)
    return y + height/2

def draw_G(width, height,x,y, ax, color="orange",lw=2.0,alpha=0.5):
    x += width/2
    y += height/2
    ts = np.linspace(3.85*math.pi/2., math.pi/6., 100)
    xs = [math.cos(t)*width/2+x for t in ts]
    ys = [math.sin(t)*height/2+y for t in ts]
    ax.plot(xs,ys,lw=lw, color=color,alpha=alpha)
    ptl = x, y-height/2
    ptu = xs[0],ys[0]
    ax.plot([ptu[0],ptu[0]], [ptu[1], ptl[1]],lw=lw, color=color,alpha=alpha)
    ax.plot([ptu[0]-width*0.75/2,ptu[0]], [ptu[1], ptu[1]],lw=lw, color=color,alpha=alpha)
    return y + height/2

def draw_seq_log(PSSM,title, AX=None,lw=1.0,alpha=1):
    ax = AX
    if AX is None:
        F = plt.figure(facecolor="white")
        ax = plt.gca()
    ax.set_title(title, fontsize=10)
    A = np.array(PSSM)
    for i in range(A.shape[0]):
        A[i,:] /= sum(A[i,:])
        methods = (draw_A, draw_U, draw_G, draw_C)
        x = i
        width = 0.9
        y = 0.0
        I = 2+sum([A[i,l]*math.log(A[i,l],2)  for l in range(4)])
        for j in range(4):
            d = methods[j]
            h = I*A[i,j]
            #print(h)
            y = d(width, h, x, y, ax, lw=lw, alpha=alpha)
    ax.set_xlim(0,A.shape[0])
    ax.set_xticks([])
    ax.set_xlabel("Position",fontsize=10)
    ax.set_ylabel("bits",fontsize=10)
    despine(ax,bottom=True)
    ax.yaxis.grid(True)
    ax.set_ylim(0,2)
    if AX is None:
        plt.show()

class logo:
    def __init__(self, matrix, name=""):
        self.m      = matrix
        self.name   = name
    def draw(self ,ax=None):
        draw_seq_log(self.m,self.name, AX=ax)

def main():
    A=[[ 0.356611,0.085184,0.455087,0.103118],
       [ 0.217352,0.033832,0.616781,0.132036],
       [ 0.411544,0.053398,0.486295,0.048762],
       [ 0.026776,0.902577,0.020027,0.050619],
       [ 0.773107,0.016504,0.055813,0.154576],
       [ 0.121956,0.027309,0.046496,0.804239],
       [ 0.010707,0.015609,0.947804,0.025881],
       [ 0.031773,0.310079,0.022036,0.636112],
       [ 0.055142,0.622140,0.063964,0.258755],
       [ 0.176130,0.035488,0.131170,0.657211],
       [ 0.095022,0.092511,0.793152,0.019314],
       [ 0.190134,0.058570,0.593589,0.157707],
       [ 0.456692,0.015034,0.473979,0.054296],
       [ 0.008153,0.976192,0.001984,0.013671],
       [ 0.767356,0.030153,0.037623,0.164869],
       [ 0.170081,0.032247,0.065448,0.732224],
       [ 0.036761,0.017177,0.926087,0.019975],
       [ 0.059032,0.416386,0.015859,0.508723],
       [ 0.114088,0.584113,0.060812,0.240987],
       [ 0.224960,0.355487,0.071233,0.348320]]
    #L        = logo(A, name="P53")
    #F        = plt.figure(figsize=(12,10),facecolor="white")
    #ax       = plt.gca()
    #L.draw(ax=ax)
    #plt.show()
if __name__ == "__main__":
    main()<jupyter_output><empty_output><jupyter_text>Extended project: Binding motifs & sequence logosTo visualize the patterns learned by each 1st layer filters, we generate a PPM (position probability matrix) derived from the filter's response to actual sequences[1].1. All positive sequences from training and test set are fed through the convolutional and rectifying (ReLU) operations of the model.2. Then, all the sequences that passed the ReLU threshold (activation>0 for one or more positions in activation map) are aligned.3. From the aligned sequences (assembly), a PPM is computed by counting the occurrences of bases in each position, which is then transformed to a sequence logo.As a result, we get 1 sequence motif per 1 convolutional filter.Below are sequence logos of learned filters:**Filter length 8**8-1![1](https://drive.google.com/uc?export=download&id=11_VIrpsoGi8LXwViwcdMgf7Hy0NaXOn6)8-2![2](https://drive.google.com/uc?export=download&id=1N52MgzFdRzbKx4Eee9q5osiAWf4ewXSK)**Filter length 16**16-1![3](https://drive.google.com/uc?export=download&id=1HVqxLAKE9_uJPE15ENLhvzGNZXK6oR4I)16-2![4](https://drive.google.com/uc?export=download&id=1x21aEjqTAhXar1cFfeppG3WIwbjewqIo)**Filter length 32**32-1![5](https://drive.google.com/uc?export=download&id=1CEuQhbOtgghPifaBNDOUojuCjPDDJYK7)32-2![6](https://drive.google.com/uc?export=download&id=1NjuAwQbJy4Jfx1f8psVQApu3mDONAttN)Let's compare the sequences to RBP-binding motifs of ELAVL1 known in literature.**Known motif**CISBP-RNA database![literature](https://drive.google.com/uc?export=download&id=1d0ckS1s8Tu9Nowgo39NXCDpygQPUwoBA)HITS-CLIP database![literature2](https://drive.google.com/uc?export=download&id=1eWzhflUmpcPfSihwTvaM-ghwLHyLqpdL)As seen from sequences above, the CNN filters learn to detect short sequences with frequent occurence of U / AU neucleotides (8-1, 8-2, 16-1, 16-2, 32-1), which is in fact similar to RBP-binding motif of ELAVL1 known from experiment.This implies that filters of CNN optimized for RBP-binding prediction is able to capture biologically important sequence motifs in RBP-binding process.**Biological meaning of this binding motif is discussed later in project2_ALKBH5.ipynb!**<jupyter_code># gather positive RBP sequences
positive = np.reshape(np.asarray(np.where(test_label[:,1])),-1) # use only positive test cases
bind_input = test_input[positive,:,:]
positive = np.reshape(np.asarray(np.where(train_label[:,1])),-1) # use only positive train cases
bind_input = np.concatenate((bind_input,train_input[positive,:,:]),axis=0)

# function for PPM computation and sequence logo generation
def gen_logo(length, w, b, logo_width):
  # generating sequence logos
  # n: filter length
  # w: filters as weight matrices
  # b: filter biases
  # logo_width: sequence logo width
  n = 2 # number of filters per each length
  n_test = bind_input.shape[0] # number of test data
  thr = 0
  #thr = np.sum(np.max(w,1),0)*0.7

  w_2d = np.reshape(w,(-1,n))
  seq = np.zeros((n_test,length,4,n)) # subsequences yielding max response for a filter
  keep = np.zeros((n_test,n)) # sequences to keep / discard

  # Align sequences for each filters
  for ii in range(n_test):
    img = bind_input[ii,:,:] # an input image, 200*4
    s = np.zeros((200-length,n)) # sliding activation map, slide*n
    for jj in range(200-length): # sliding convolution, 'valid'
      img_1d = np.reshape(img[jj:jj+length,:],(1,-1)) # 1*24
      s[jj,:] = np.matmul(img_1d,w_2d) # convolution operation
      s[jj,:] += b # bias addition
    
    max_s = np.max(s,0) # compute maximum activation of each filters
    keep[ii,:] = max_s>thr # determine exclusion from assembly
    # select max scoring subsequence for each filter
    max_t = np.argmax(s,0)
    for ff in range(n): seq[ii,:,:,ff] = img[max_t[ff]:max_t[ff]+length,:]

  # Select positive-responding sequences and create PPM of each filter
  for ff in range(n):
    # compute PPM of each filter
    tmp = np.reshape(np.asarray(np.where(keep[:,ff])),-1) # sequences to exclude / include
    ppm = np.mean(seq[tmp,:,:,ff],0) # compute PPM from aligned assembly
    ppm[np.where(ppm==0)] += 1e-10 # for numerical stability

    # generate sequence logo
    L = logo(ppm, name=str(ff+1))
    F = plt.figure(figsize=(logo_width,3),facecolor="white")
    ax = plt.gca()
    L.draw(ax=ax)
    plt.show()
gen_logo(length=8, w=w1, b=b1, logo_width = 5)
gen_logo(length=16, w=w2, b=b2, logo_width = 10)
gen_logo(length=32, w=w3, b=b3, logo_width = 15)<jupyter_output><empty_output>
{'fraction_non_alphanumeric': 0.13903811716551548, 'fraction_numerical': 0.05222989681930588, 'mean_word_length': 4.675054439874184, 'pattern_counts': {'":': 1, '<': 31, '<?xml version=': 0, '>': 43, 'https://': 15, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 13, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '7', 'dir': 'jupyter-structured-clean-dedup', 'id': '471459', 'license': 'MIT', 'n_tokens_mistral': 8981, 'n_tokens_neox': 8260, 'n_words': 2139, 'path': 'Project#2/project2_ELAVL1.ipynb', 'repo_name': 'bugkjw/bis438-bioinformatics', 'size': '23466'}