text:
#-------------------------------------------------------------------------
# Copyright (c) Microsoft Corporation.  All rights reserved.
# Licensed under the MIT License.
#--------------------------------------------------------------------------

import logging
import onnx
import sys
import argparse
import numpy as np
from collections import deque
from onnx import ModelProto, TensorProto, numpy_helper, helper
from onnx_model_bert import BertOnnxModel

logger = logging.getLogger(__name__)

class BertOnnxModelTF(BertOnnxModel):
    def __init__(self, model, num_heads, hidden_size):
        super().__init__(model, num_heads, hidden_size)

    def remove_identity(self):
        nodes_to_remove = []
        for node in self.nodes():
            if node.op_type == 'Identity':
                if not self.find_graph_output(node.output[0]):
                    self.replace_input_of_all_nodes(node.output[0], node.input[0])
                    nodes_to_remove.append(node)
        self.remove_nodes(nodes_to_remove)
        logger.info(f"Removed Identity count: {len(nodes_to_remove)}")

    def match_mask_path(self, add_or_sub_before_softmax):
        mask_nodes = self.match_parent_path(add_or_sub_before_softmax, ['Mul', 'Sub', 'Reshape', 'Cast'],
                                            [1, None, 1, 0])
        if mask_nodes is not None:
            return mask_nodes

        mask_nodes = self.match_parent_path(add_or_sub_before_softmax, ['Mul', 'Sub', 'Cast', 'Slice', 'Unsqueeze'],
                                            [1, 0, 1, 0, 0])
        if mask_nodes is not None:
            return mask_nodes

        mask_nodes = self.match_parent_path(add_or_sub_before_softmax, ['Mul', 'Sub', 'Cast', 'Unsqueeze', 'Unsqueeze'],
                                            [1, None, 1, 0, 0])

        return mask_nodes

    def get_2d_initializers_from_parent_subgraphs(self, current_node):
        """
        Find initializers that is 2D. Returns a dictionary with name as key and shape as value.
        """
        parent_nodes = self.get_parent_subgraph_nodes(current_node, [])
        initializers = {}
        for node in parent_nodes:
            for input in node.input:
                initializer = self.get_initializer(input)
                if initializer:
                    temp = numpy_helper.to_array(initializer)
                    if len(temp.shape) == 2:
                        initializers[initializer.name] = temp.shape

        return initializers

    def find_segment_ids(self, segment_embedding, input_ids):
        input_name_to_nodes = self.input_name_to_nodes()
        if segment_embedding not in input_name_to_nodes:
            return None

        nodes = input_name_to_nodes[segment_embedding]
        if len(nodes) != 1:
            return None

        graph_inputs = self.get_graph_inputs(nodes[0], recursive=True)
        if len(graph_inputs) > 1:
            print("Found multiple candidates of segment_ids", graph_inputs)
            return None
        # Find segment ids in graph inputs. The segment id input must not be the same as input_ids.
        if len(graph_inputs) == 1 and graph_inputs[0] != input_ids:
            return graph_inputs[0]

        # If the segment id candidate is the same as the input_ids, try to assign alternative segment ids and simplify the graph if needed.
        segment_ids = nodes[0].input[1]
        _, segment_id_path, _ = self.match_parent_paths(
            nodes[0], [(["ConstantOfShape", "Cast", "Concat", "Slice", "Cast", "Shape"], [1, 0, 0, 0, 0, 0]),
                       (["ConstantOfShape", "Cast", "Concat", "Unsqueeze", "Squeeze", "Slice", "Cast", "Shape"
                         ], [1, 0, 0, 0, 0, 0, 0, 0])], None)

        if segment_id_path and input_ids and input_ids == segment_id_path[-1].input[0]:
            logger.debug("Simplify semgent id path...")
            constantofshape_node = segment_id_path[0]
            graph_name = self.get_graph_by_node(constantofshape_node).name
            self.add_node(helper.make_node('Shape', inputs=[input_ids], outputs=["input_shape"]), graph_name)
            constantofshape_value = helper.get_attribute_value(constantofshape_node.attribute[0])
            self.add_node(
                helper.make_node('ConstantOfShape',
                                 inputs=["input_shape"],
                                 outputs=["zeros_for_input_shape"],
                                 value=constantofshape_value), graph_name)
            segment_ids = "zeros_for_input_shape"
        return segment_ids

    def find_input_ids(self, word_embedding):
        input_name_to_nodes = self.input_name_to_nodes()
        if word_embedding not in input_name_to_nodes:
            return None

        nodes = input_name_to_nodes[word_embedding]
        if len(nodes) != 1:
            return None

        graph_inputs = self.get_graph_inputs(nodes[0], recursive=True)
        if len(graph_inputs) == 1:
            return graph_inputs[0]

        print("Found multiple candidates of input_ids", graph_inputs)
        return None

    def find_mask_input(self, excluded_graph_inputs):
        for node in self.nodes():
            if node.op_type == 'Softmax':
                mask_path = self.match_parent_path(node, ['Add', 'Mul', 'Sub', 'Cast', 'Slice', 'Unsqueeze'],
                                                   [0, 1, None, 1, 0, 0])
                if mask_path is None:
                    continue
                add_node, mul_node, sub_node, cast_node, slice_node, unsqueeze_node = mask_path
                if self.has_constant_input(mul_node, -10000) and self.has_constant_input(sub_node, 1):
                    graph_inputs = self.get_graph_inputs(sub_node, recursive=True)
                    inputs = [input for input in graph_inputs if input not in excluded_graph_inputs]
                    if len(inputs) > 1:
                        print("Found multiple candidates of mask input", inputs)
                        return None
                    if len(inputs) == 1:
                        return inputs[0]
                    # Duplicated input found. Try to simplify the graph.
                    path_to_be_simplified = self.match_parent_path(
                        mask_path[-1],
                        ["ConstantOfShape", "Cast", "Concat", "Unsqueeze", "Squeeze", "Slice", "Cast", "Shape"],
                        [0, 0, 0, 0, 0, 0, 0, 0])
                    duplicated_inputs = [input for input in graph_inputs if input in excluded_graph_inputs]
                    # Simplify graph for dynamic axes.
                    if path_to_be_simplified and duplicated_inputs and len(
                            duplicated_inputs) == 1 and duplicated_inputs[0] == path_to_be_simplified[-1].input[0]:
                        logger.debug("Simplify semgent id path...")
                        constantofshape_node = path_to_be_simplified[0]
                        constantofshape_value = helper.get_attribute_value(constantofshape_node.attribute[0])
                        graph_name = self.get_graph_by_node(constantofshape_node).name
                        self.add_node(
                            helper.make_node('Shape', inputs=[duplicated_inputs[0]], outputs=["input_shape_for_mask"]),
                            graph_name)
                        self.add_node(
                            helper.make_node('ConstantOfShape',
                                             inputs=["input_shape_for_mask"],
                                             outputs=[unsqueeze_node.input[0]],
                                             value=constantofshape_value), graph_name)
                    return unsqueeze_node.input[0]
        return None

    def create_embedding_subgraph(self, normalize_node, word_embedding, segment_embedding, position_embedding):
        input_ids = self.find_input_ids(word_embedding)
        if input_ids is None:
            logger.info("Failed to find input_ids. Cannot fuse embedding layer.")
            return False

        segment_ids = self.find_segment_ids(segment_embedding, input_ids)
        if segment_ids is None:
            logger.info("Failed to find segment_ids. Cannot fuse embedding layer.")
            return False

        mask_input = self.find_mask_input([segment_ids, input_ids])
        if mask_input is None:
            logger.info("Failed to find input_mask. Cannot fuse embedding layer.")
            return False

        self.bert_inputs = [input_ids, segment_ids, mask_input]

        mask_index = self.create_node_name('mask_index')
        self.attention_mask.set_mask_indice(mask_input, mask_index)

        if self.find_graph_input(input_ids).type.tensor_type.elem_type != TensorProto.INT32:
            casted, input_ids = self.utils.cast_graph_input_to_int32(input_ids)

        if self.find_graph_input(segment_ids):
            casted, segment_ids = self.utils.cast_graph_input_to_int32(segment_ids)
        else:
            segment_ids, segment_id_cast_node = self.utils.cast_input_to_int32(segment_ids)

        if self.find_graph_input(mask_input):
            casted, mask_input = self.utils.cast_graph_input_to_int32(mask_input)
        else:
            mask_input, mask_input_cast_node = self.utils.cast_input_to_int32(mask_input)

        embed_output = self.create_node_name('embed_output')
        embed_node = onnx.helper.make_node(
            'EmbedLayerNormalization',
            inputs=[
                input_ids,
                segment_ids,
                word_embedding,
                position_embedding,
                segment_embedding,
                normalize_node.input[1],  # gamma
                normalize_node.input[2],  # beta
                mask_input
            ],
            outputs=[embed_output, mask_index],
            name="EmbedLayer")
        embed_node.domain = "com.microsoft"
        self.replace_input_of_all_nodes(normalize_node.output[0], embed_output)
        self.add_node(embed_node, self.get_graph_by_node(normalize_node).name)

    def process_embedding(self):
        """
        Automatically detect word, segment and position embeddings.
        """
        logger.info("start processing embedding layer...")
        output_name_to_node = self.output_name_to_node()

        layer_norm_nodes = self.get_nodes_by_op_type("LayerNormalization")
        for layer_norm_node in layer_norm_nodes:
            pos_embed_path = self.match_parent_path(layer_norm_node, ['Add', 'Reshape', 'Slice'], [0, 1, 0],
                                                    output_name_to_node)
            if pos_embed_path is None:
                continue

            add_node, reshape_node, slice_node = pos_embed_path
            initializer = self.get_initializer(slice_node.input[0])
            if initializer is None:
                continue

            temp = numpy_helper.to_array(initializer)
            if len(temp.shape) == 2:
                logger.info("Found position embedding. name:{}, shape:{}".format(initializer.name, temp.shape))
                position_embedding = initializer.name
            else:
                logger.info("Failed to find position embedding. name:{}, shape:{}".format(initializer.name, temp.shape))
                return

            first_parent = self.get_parent(add_node, 0, output_name_to_node)
            if first_parent is not None and first_parent.op_type == "Add":
                embeddings = self.get_2d_initializers_from_parent_subgraphs(first_parent)
                if len(embeddings) != 2:
                    logger.warning(
                        "Failed to find two embeddings (word and segment) from Add node. Found {}".format(embeddings))
                    return

                word_embedding = None
                segment_embedding = None
                for name, shape in embeddings.items():
                    if shape[0] == 2:
                        segment_embedding = name
                        logger.info("Found segment embedding. name:{}, shape:{}".format(name, shape))
                    else:
                        word_embedding = name
                        logger.info("Found words embedding. name:{}, shape:{}".format(name, shape))

                if word_embedding is None or segment_embedding is None:
                    logger.info("Failed to find both word and segment embedding")
                    return

                logger.info("Create Embedding node")
                self.create_embedding_subgraph(layer_norm_node, word_embedding, segment_embedding, position_embedding)
                # Prune graph to remove those original embedding nodes.
                self.prune_graph()
                break

    def check_attention_input(self, matmul_q, matmul_k, matmul_v, parent, output_name_to_node):
        for x in [matmul_q, matmul_k, matmul_v]:
            root_input = x.input[0]
            root_node = output_name_to_node[root_input]
            if root_node == parent:
                continue
            logger.debug(f"Check attention input failed:{root_input}, {parent.output[0]}")
            return False

        return True

    def fuse_attention(self):
        output_name_to_node = self.output_name_to_node()

        nodes_to_remove = []
        attention_count = 0

        start_nodes = []
        skip_layer_norm_nodes = self.get_nodes_by_op_type("SkipLayerNormalization")
        layer_norm_nodes = self.get_nodes_by_op_type("LayerNormalization")
        # Sometimes we can not fuse skiplayernormalization since the add before layernorm has an output that used by nodes outside skiplayernorm
        # Conceptually we treat add before layernorm as skiplayernorm node since they share the same pattern
        start_nodes.extend(skip_layer_norm_nodes)
        start_nodes.extend(layer_norm_nodes)

        for normalize_node in start_nodes:
            graph_name = self.get_graph_by_node(normalize_node).name
            # SkipLayerNormalization has two inputs, and one of them is the root input for attention.
            if normalize_node.op_type == 'LayerNormalization':
                add_before_layernorm = self.match_parent(normalize_node, 'Add', 0)
                if add_before_layernorm is not None:
                    normalize_node = add_before_layernorm
                else:
                    continue
            parent = self.get_parent(normalize_node, 1)
            if parent is None or parent.op_type not in ["SkipLayerNormalization", "LayerNormalization", "Reshape"]:
                parent = self.get_parent(normalize_node, 0)
                if parent is None or parent.op_type not in ["SkipLayerNormalization", "LayerNormalization", "Reshape"]:
                    logger.debug("Failed to match parent of normalize_node")
                    continue

            qkv_nodes = self.match_parent_path(normalize_node, ['Add', 'MatMul', 'Reshape', 'Transpose', 'MatMul'],
                                               [0, 0, 0, 0, 0])
            if qkv_nodes is None:
                qkv_nodes = self.match_parent_path(normalize_node, ['MatMul', 'Reshape', 'Transpose', 'MatMul'],
                                                   [1, 0, 0, 0])
                if qkv_nodes is None:
                    qkv_nodes = self.match_parent_path(normalize_node, ['Add', 'Einsum', 'Einsum'], [0, 0, 0])
                    if qkv_nodes is None:
                        logger.debug("Failed to match qkv nodes")
                        continue

            matmul_qkv = qkv_nodes[-1]
            v_nodes = self.match_parent_path(matmul_qkv, ['Transpose', 'Reshape', 'Add', 'MatMul'], [1, 0, 0, 0])
            if v_nodes is None:
                v_nodes = self.match_parent_path(matmul_qkv, ['Add', 'Einsum'], [1, 0])
                if v_nodes is None:
                    logger.debug("Failed to match v path")
                    continue

            add_v = v_nodes[-2]
            matmul_v = v_nodes[-1]
            qk_nodes = self.match_parent_path(matmul_qkv, ['Softmax', 'Add', "Mul", 'MatMul'], [0, 0, 0, 0])
            if qk_nodes is None:
                qk_nodes = self.match_parent_path(matmul_qkv, ['Softmax', 'Add', 'Einsum'], [0, 0, 0])
                if qk_nodes is None:
                    logger.debug("Failed to match qk_paths")
                    continue
            matmul_qk = qk_nodes[-1]

            q_nodes = self.match_parent_path(matmul_qk, ['Transpose', 'Reshape', 'Add', 'MatMul'], [0, 0, 0, 0])
            if q_nodes is None:
                q_nodes = self.match_parent_path(matmul_qk, ['Add', 'Einsum'], [0, 0])
                if q_nodes is None:
                    logger.debug("Failed to match q path")
                    continue

            add_q = q_nodes[-2]
            matmul_q = q_nodes[-1]

            k_nodes = self.match_parent_path(matmul_qk, ['Transpose', 'Reshape', 'Add', 'MatMul'], [1, 0, 0, 0])
            if k_nodes is None:
                k_nodes = self.match_parent_path(matmul_qk, ['Mul', 'Add', 'Einsum'], [1, 0, 0])
                if k_nodes is None:
                    logger.debug("Failed to match k path")
                    continue
            add_k = k_nodes[-2]
            matmul_k = k_nodes[-1]

            mask_nodes = self.match_mask_path(qk_nodes[1])

            if mask_nodes is None:
                logger.debug("Cannot find mask_nodes.")
                continue

            if not self.has_constant_input(mask_nodes[1], 1):
                logger.debug("Sub node expected to have an input with constant value 1.0.")
                continue

            # add a squeeze node to convert a 3-d mask to 2-d
            squeeze_node = self.match_parent_path(mask_nodes[-1], ['Squeeze'], [0]) or self.match_parent_path(
                mask_nodes[-1], ['Expand'], [0])
            squeeze_node_name = "Squeeze_3d_to_2d_mask"
            squeeze_output_name = squeeze_node_name + "_output"
            if squeeze_node is None and len(mask_nodes) == 5 and self.find_graph_input(mask_nodes[-1].input[0]) is None:
                mask_input = mask_nodes[-1].input[1]
                self.add_node(
                    helper.make_node("Squeeze", [mask_input], [squeeze_output_name], squeeze_node_name, axes=[1]),
                    graph_name)
                mask_nodes[-1].input[0] = squeeze_output_name

            is_same_root = self.check_attention_input(matmul_q, matmul_k, matmul_v, parent, output_name_to_node)
            if is_same_root:
                mask_index = self.attention_mask.process_mask(mask_nodes[-1].input[0])
                logger.debug("Create an Attention node.")

                # For tf models, q and v are flipped.
                attention_node = self.attention_fusion.create_attention_node(mask_index, matmul_k, matmul_q, matmul_v,
                                                                             add_k, add_q, add_v, self.num_heads,
                                                                             self.hidden_size, parent.output[0],
                                                                             qkv_nodes[2].output[0], None)
                if attention_node is None:
                    continue

                if qkv_nodes[1].op_type == 'Einsum':
                    # add reshape before einsum
                    tensor = helper.make_tensor(name=qkv_nodes[1].name + "_newshape",
                                                data_type=TensorProto.INT64,
                                                dims=[4],
                                                vals=np.int64(
                                                    [[0, 0, self.num_heads,
                                                      int(self.hidden_size / self.num_heads)]]).tobytes(),
                                                raw=True)
                    self.add_initializer(tensor, graph_name)
                    reshape_ = helper.make_node("Reshape",
                                                inputs=[attention_node.output[0], qkv_nodes[1].name + "_newshape"],
                                                outputs=[qkv_nodes[1].name + "_reshape_output"],
                                                name=qkv_nodes[1].name + "_reshape")
                    qkv_nodes[1].input[0] = qkv_nodes[1].name + "_reshape_output"
                    self.add_node(reshape_, graph_name)
                if parent.op_type == 'Reshape':
                    # Temporary work around: we require the skiplayernorm and attention op be fed with 3-d input
                    hidden_size = numpy_helper.to_array(self.get_initializer(parent.input[1]))[1]
                    tensor = helper.make_tensor(name=parent.name + "_modified",
                                                data_type=TensorProto.INT64,
                                                dims=[3],
                                                vals=np.int64([[1, -1, hidden_size]]).tobytes(),
                                                raw=True)
                    self.add_initializer(tensor, graph_name)
                    parent.input[1] = parent.name + "_modified"

                self.add_node(attention_node, graph_name)
                attention_count += 1

                nodes_to_remove.extend(qkv_nodes[2:])
                nodes_to_remove.extend(qk_nodes)
                nodes_to_remove.extend(q_nodes)
                nodes_to_remove.extend(k_nodes)
                nodes_to_remove.extend(v_nodes)
                nodes_to_remove.extend(mask_nodes)
            else:
                logger.debug("Root node not matched.")
                continue
        self.remove_nodes(nodes_to_remove)
        self.update_graph()
        logger.info(f"Fused Attention count:{attention_count}")

    def preprocess(self):
        self.remove_identity()
        self.process_embedding()
        self.skip_reshape()

    def skip_reshape(self):
        count = 0
        reshape_nodes = self.get_nodes_by_op_type("Reshape")
        for reshape_node in reshape_nodes:
            parent = self.get_parent(reshape_node, 0)
            if parent is not None and parent.op_type == "Reshape":
                reshape_node.input[0] = parent.input[0]
                count += 1

        if count > 0:
            logger.info(f"Skip consequent Reshape count: {count}")

    def remove_reshape_before_first_attention(self):
        attention_nodes = self.get_nodes_by_op_type("Attention")
        for attention_node in attention_nodes:
            path = self.match_parent_path(attention_node, ['Reshape', 'EmbedLayerNormalization'], [0, 0])
            if path is None:
                continue
            logger.info("Remove Reshape before first Attention node.")
            reshape, _ = path
            self.replace_input_of_all_nodes(reshape.output[0], reshape.input[0])
            self.remove_node(reshape)
            break

    def postprocess(self):
        self.remove_reshape_before_first_attention()
        self.prune_graph()

max_stars_repo_name: raikonenfnu/transformer-benchmarks

max_issues_repo_name: raikonenfnu/transformer-benchmarks

max_forks_repo_name: raikonenfnu/transformer-benchmarks


text:
from __future__ import print_function, division
import visdom
import torch
import torch.nn as nn
from torch.autograd import Variable
import numpy as np
import torchvision
from torchvision import datasets, models, transforms
import torch.utils.model_zoo as model_zoo
import time
from os.path import join, exists, isdir
from os import makedirs, listdir
import argparse
import copy
from PIL import Image
from matplotlib import pyplot as plt
from collections import namedtuple
######################################################################################################################################
__all__ = [
    'vgg16_bn', 'vgg16'
]

model_urls = {
    'vgg16': 'https://download.pytorch.org/models/vgg16-397923af.pth',
    'vgg16_bn': 'https://download.pytorch.org/models/vgg16_bn-6c64b313.pth'
}

class VGG(nn.Module):

    def __init__(self, features, num_classes=1000, init_weights=True, dropout=[0.5,0.5]):
        super(VGG, self).__init__()
        self.features = features
        self.avgpool = nn.AdaptiveAvgPool2d((7, 7))
        self.classifier = nn.Sequential(
            nn.Linear(512 * 7 * 7, 4096),
            nn.ReLU(True),
            nn.Dropout(dropout[0]),
            nn.Linear(4096, 4096),
            nn.ReLU(True),
            nn.Dropout(dropout[1]),
            nn.Linear(4096, num_classes),
        )
        if init_weights:
            self._initialize_weights()

    def forward(self, x):
        x = self.features(x)
        x = self.avgpool(x)
        x = x.view(x.size(0), -1)
        x = self.classifier(x)
        return x

    def _initialize_weights(self):
        for m in self.modules():
            if isinstance(m, nn.Conv2d):
                nn.init.kaiming_normal_(m.weight, mode='fan_out', nonlinearity='relu')
                if m.bias is not None:
                    nn.init.constant_(m.bias, 0)
            elif isinstance(m, nn.BatchNorm2d):
                nn.init.constant_(m.weight, 1)
                nn.init.constant_(m.bias, 0)
            elif isinstance(m, nn.Linear):
                nn.init.normal_(m.weight, 0, 0.01)
                nn.init.constant_(m.bias, 0)

def make_layers(cfg, batch_norm=False):
    layers = []
    in_channels = 3
    for v in cfg:
        if v == 'M':
            layers += [nn.MaxPool2d(kernel_size=2, stride=2)]
        else:
            conv2d = nn.Conv2d(in_channels, v, kernel_size=3, padding=1)
            if batch_norm:
                layers += [conv2d, nn.BatchNorm2d(v), nn.ReLU(inplace=True)]
            else:
                layers += [conv2d, nn.ReLU(inplace=True)]
            in_channels = v
    return nn.Sequential(*layers)

cfg = {
    'A': [64, 'M', 128, 'M', 256, 256, 'M', 512, 512, 'M', 512, 512, 'M'],
    'B': [64, 64, 'M', 128, 128, 'M', 256, 256, 'M', 512, 512, 'M', 512, 512, 'M'],
    'D': [64, 64, 'M', 128, 128, 'M', 256, 256, 256, 'M', 512, 512, 512, 'M', 512, 512, 512, 'M'],
    'E': [64, 64, 'M', 128, 128, 'M', 256, 256, 256, 256, 'M', 512, 512, 512, 512, 'M', 512, 512, 512, 512, 'M'],
}

def vgg16(pretrained=False, **kwargs):
    """VGG 16-layer model (configuration "D")

    Args:
        pretrained (bool): If True, returns a model pre-trained on ImageNet
    """
    if pretrained:
        kwargs['init_weights'] = False
    model = VGG(make_layers(cfg['D']), **kwargs)
    if pretrained:
        model.load_state_dict(model_zoo.load_url(model_urls['vgg16']))
    return model

def vgg16_bn(pretrained=False, **kwargs):
    """VGG 16-layer model (configuration "D") with batch normalization
    Args:
        pretrained (bool): If True, returns a model pre-trained on ImageNet
    """
    if pretrained:
        kwargs['init_weights'] = False
    model = VGG(make_layers(cfg['D'], batch_norm=True), **kwargs)
    if pretrained:
        model.load_state_dict(model_zoo.load_url(model_urls['vgg16_bn']))
    return model

def plot_kernels(tensor, num_cols=8):
    num_kernels = tensor.shape[0]
    num_rows = 1 + num_kernels // num_cols
    fig = plt.figure(figsize=(num_cols,num_rows))
    for i in range(tensor.shape[0]):
        ax1 = fig.add_subplot(num_rows,num_cols,i+1)
        ax1.imshow(tensor[i], cmap='gray')
        ax1.axis('off')
        ax1.set_xticklabels([])
        ax1.set_yticklabels([])

    plt.subplots_adjust(wspace=0.1, hspace=0.1)
    return plt

######################################################################################################################################
######################################################################################################################################

parser = argparse.ArgumentParser(description='Visualize weights of VGG model using Visdom.')
parser.add_argument('model_name', metavar='model-name', help='Names of dataset on which model was trained. It will be used as suffix for model directory.', type=str)
parser.add_argument('model_dir', metavar='model-dir', help='Directory where model for dataset is to be found.', type=str)
parser.add_argument('--batchnorm', help='Whether Batch Normalization was used.', action='store_true')
parser.add_argument('--metrics', help='Whether Metrics should be displayed.', action='store_true')
parser.add_argument('--visualize', help='Whether to visualize', action='store_true')
parser.add_argument('-l', '--layer', help='Relu Layer number. Default 0. choices are 0-4.', nargs='?', type=int, const=0, default=0, choices=range(5))
parser.add_argument('-f', '--filter', help='Filter number. Default 0.', nargs='?', type=int, const=0, default=0)
args = parser.parse_args()

if not isdir(join(args.model_dir, args.model_name)):
    raise Exception('Given model directory does not exist')

model_dir = join(args.model_dir, args.model_name)

if args.metrics:
    vis = visdom.Visdom(server='http://10.4.16.22', env="metrics")
    metrics = {'tr': dict(), 'vl': dict()}
    names = set()
    print(listdir(model_dir))
    files = [x for x in listdir(model_dir) if x.startswith('vgg-tr-') or x.startswith('vgg-vl-')]
    print(files)
    for metric in files:
        name = metric.split('-')[-1].split('.')[0]
        phase = metric.split('-')[-2]
        names.add(name)
        file = join(model_dir, metric)
        y = [float(value if not value.startswith('tensor') else value.split('(')[1].split(',')[0]) for value in open(file)]
        x = [x for x in range(len(y))]
        metrics[phase][name] = {'x': x, 'y': y}
    for name in names:
        for phase in metrics.keys():
            if name in metrics[phase]:
                plt.plot(metrics[phase][name]['x'], metrics[phase][name]['y'], label = phase)
        plt.xlabel('Epochs')
        plt.ylabel(name.upper())
        plt.title('Model ' + name.upper())
        plt.legend()
        vis.matplot(plt)
        plt.close()

######################################################################################################################################
# GPU and Creating Model
if args.visualize:
    class Vgg16(torch.nn.Module):
        def __init__(self, model):
            super(Vgg16, self).__init__()
            features = list(model.features)[:30]
            self.features = nn.ModuleList(features).eval()

        def forward(self, x):
            results = []
            for ii, model in enumerate(self.features):
                x = model(x)
                if ii in {3, 8, 15, 22, 29}:
                    results.append(x)

            vgg_outputs = namedtuple("VggOutputs", ['relu1_2', 'relu2_2', 'relu3_3', 'relu4_3', 'relu5_3'])
            return vgg_outputs(*results)

    device = torch.device("cuda:0" if torch.cuda.is_available() else "cpu")

    model_conv = vgg16_bn(pretrained=True, dropout=[0.3,0.3]) if args.batchnorm else vgg16(pretrained=True, dropout=[0.3,0.3])
    new_model_conv = (list(model_conv.classifier.children())[:-1])
    new_model_conv.append(nn.Linear(4096, 2))
    model_conv.classifier = nn.Sequential(*new_model_conv)
    model_conv = nn.DataParallel(model_conv)
    model_conv.load_state_dict(torch.load(join(model_dir, 'vgg_wei_err.pth'), device))
    model = Vgg16(model_conv.module.double())
    print(model)

    vis = visdom.Visdom(server='http://10.4.16.22')

    img = None
    relus = dict()
    for x in ['cancer', 'normal']:
    	img = Image.open('../input_images/' + ('lung' if args.model_name.startswith('lung') else 'prostate') + '_' + x + '_01.jpg').resize((224,224))
    	vis.image(np.rollaxis(np.array(img), 2, 0), opts={'title': x})
    	inp = np.ndarray((1, 3, 224, 224), dtype=np.double)
    	inp[0] = np.rollaxis(np.array(img), 2, 0)
    	relus[x] = model(torch.from_numpy(inp))

    plot = dict()
    # Properties window
    properties = [
        {'type': 'select', 'name': 'Block', 'value': 0, 'values': [0, 1, 2, 3, 4]},
        {'type': 'select', 'name': 'Filter', 'value': 0, 'values': [0, 1, 2, 3, 4]}
    ]

    properties_window = vis.properties(properties, opts={'title': 'Properties'})

    def properties_callback(event):
    	if event['event_type'] == 'PropertyUpdate':
    		prop_id = event['propertyId']
    		value = event['value']
    		properties[prop_id]['value'] = int(value)
    		vis.properties(properties, win=properties_window)
    		visualize(properties[0]['value'], properties[1]['value'])

    vis.register_event_handler(properties_callback, properties_window)

    def visualize(block_num, filter_num):
        layer = list([2, 7, 14, 21, 28])[block_num]
        print('found layer')
        trained_feature = np.uint8((list(model_conv.module.features.children())[layer]).weight.double().data.numpy()*255)
        print(trained_feature.shape)
        properties[1]['values'] = [x for x in range(trained_feature.shape[0])]
        vis.properties(properties, win=properties_window)
        plot['weights'] = vis.matplot(plot_kernels(trained_feature[filter_num] if layer else trained_feature, int(np.sqrt(trained_feature.shape[0]))), win = plot['weights'] if "weights" in plot else None, opts={'title': 'Filter Components'})
        print('plotted weights')
        for x in ['cancer', 'normal']:
            block = relus[x][block_num].detach().numpy()
            print(block.shape)
            imgs = block[0]
            plot[x] = vis.matplot(plot_kernels(imgs, int(np.sqrt(imgs.shape[0]))), win=plot[x] if x in plot else None, opts={'title': 'All Filter Results for ' + x})
            plot[x+"_filter"] = vis.image(np.array(Image.fromarray(np.uint8(block[0][filter_num])).resize((224,224))), win=plot[x+"_filter"] if x+"_filter" in plot else None, opts={'title': 'Specific Filter Result for ' + x})
            print('plotted for', x, 'image')

    visualize(properties[0]['value'], properties[1]['value'])

    try:
    	input = raw_input  # for Python 2 compatibility
    except NameError:
    	pass
    	input('Waiting for callbacks, press enter to quit.')

ext: py

reason: import numpy


text:
"""Fill NA."""
# import numpy as np
import pandas as pd

def fill_zero(df):
    """Fill NaN with zero."""
    df = df.fillna(0)
    return df

def fill_forward(df):
    """Fill NaN with previous values."""
    df = df.fillna(method='ffill')
    df = df.fillna(method='bfill').fillna(0)
    return df

def fill_mean(df):
    """Fill NaN with mean."""
    df = df.fillna(df.mean().fillna(0).to_dict())
    return df

def fill_median(df):
    """Fill NaN with median."""
    df = df.fillna(df.median().fillna(0).to_dict())
    return df

def rolling_mean(df, window: int = 10):
    """Fill NaN with mean of last window values."""
    df = fill_forward(df.fillna(df.rolling(window=window, min_periods=1).mean()))
    return df

def biased_ffill(df, mean_weight: float = 1):
    """Fill NaN with average of last value and mean."""
    df_mean = fill_mean(df)
    df_ffill = fill_forward(df)
    df = ((df_mean * mean_weight) + df_ffill) / (1 + mean_weight)
    return df

def fake_date_fill(df, back_method: str = 'slice'):
    """
    Return a dataframe where na values are removed and values shifted forward.

    Warnings:
        Thus, values will have incorrect timestamps!

    Args:
        back_method (str): how to deal with tails left by shifting NaN
            - 'bfill' -back fill the last value
            - 'slice' - drop any rows with any na
            - 'keepna' - keep the lagging na
    """
    df_index = df.index.to_series().copy()
    df2 = df.sort_index(ascending=False).copy()
    df2 = df2.apply(lambda x: pd.Series(x.dropna().values))
    df2 = df2.sort_index(ascending=False)
    df2.index = df_index.tail(len(df2.index))
    df2 = df2.dropna(how='all', axis=0)
    if df2.empty:
        df2 = df.fillna(0)

    if back_method == 'bfill':
        df2 = fill_forward(df2)
        return df
    elif back_method == 'slice':
        thresh = int(df.shape[1] * 0.5)
        thresh = thresh if thresh > 1 else 1
        df3 = df2.dropna(thresh=thresh, axis=0)
        if df3.empty or df3.shape[0] < 8:
            df3 = fill_forward(df2)
        else:
            df3 = fill_forward(df3)
        return df3
    elif back_method == 'keepna':
        return df2
    else:
        print('back_method not recognized in fake_date_fill')
        return df2

df_interpolate = {
    'linear': 0.1,
    'time': 0.1,
    'pad': 0.1,
    'nearest': 0.1,
    'zero': 0.1,
    'quadratic': 0.1,
    'cubic': 0.1,
    'spline': 0.1,
    'barycentric': 0.01,  # this parallelizes and is noticeably slower
    'piecewise_polynomial': 0.01,
    'spline': 0.1,
    'pchip': 0.1,
    'akima': 0.1,
    # these seem to cause more harm than good usually
    'polynomial': 0.0,
    'krogh': 0.0,
    'cubicspline': 0.0,
    'from_derivatives': 0.0,
    'slinear': 0.0,
}
df_interpolate_full = list(df_interpolate.keys())

def FillNA(df, method: str = 'ffill', window: int = 10):
    """Fill NA values using different methods.

    Args:
        method (str):
            'ffill' - fill most recent non-na value forward until another non-na value is reached
            'zero' - fill with zero. Useful for sales and other data where NA does usually mean $0.
            'mean' - fill all missing values with the series' overall average value
            'median' - fill all missing values with the series' overall median value
            'rolling mean' - fill with last n (window) values
            'ffill mean biased' - simple avg of ffill and mean
            'fake date' - shifts forward data over nan, thus values will have incorrect timestamps
            also most `method` values of pd.DataFrame.interpolate()
        window (int): length of rolling windows for filling na, for rolling methods
    """
    method = str(method).replace(" ", "_")

    if method == 'zero':
        return fill_zero(df)

    elif method == 'ffill':
        return fill_forward(df)

    elif method == 'mean':
        return fill_mean(df)

    elif method == 'median':
        return fill_median(df)

    elif method == 'rolling_mean':
        return rolling_mean(df, window=window)

    elif method == 'rolling_mean_24':
        return rolling_mean(df, window=24)

    elif method == 'ffill_mean_biased':
        return biased_ffill(df)

    elif method == 'fake_date':
        return fake_date_fill(df, back_method='slice')

    elif method in df_interpolate_full:
        df = df.interpolate(method=method, order=5).fillna(method='bfill')
        if df.isnull().values.any():
            df = fill_forward(df)
        return df

    elif method == 'IterativeImputer':
        cols = df.columns
        indx = df.index
        try:
            from sklearn.experimental import enable_iterative_imputer  # noqa
        except Exception:
            pass
        from sklearn.impute import IterativeImputer

        df = IterativeImputer(random_state=0, max_iter=100).fit_transform(df)
        if not isinstance(df, pd.DataFrame):
            df = pd.DataFrame(df)
            df.index = indx
            df.columns = cols
        return df

    elif method == 'IterativeImputerExtraTrees':
        cols = df.columns
        indx = df.index
        try:
            from sklearn.experimental import enable_iterative_imputer  # noqa
        except Exception:
            pass
        from sklearn.ensemble import ExtraTreesRegressor
        from sklearn.impute import IterativeImputer

        df = IterativeImputer(
            ExtraTreesRegressor(n_estimators=10, random_state=0),
            random_state=0,
            max_iter=100,
        ).fit_transform(df)
        if not isinstance(df, pd.DataFrame):
            df = pd.DataFrame(df)
            df.index = indx
            df.columns = cols
        return df

    elif method == 'KNNImputer':
        cols = df.columns
        indx = df.index
        from sklearn.impute import KNNImputer

        df = KNNImputer(n_neighbors=5).fit_transform(df)
        if not isinstance(df, pd.DataFrame):
            df = pd.DataFrame(df)
            df.index = indx
            df.columns = cols
        return df

    elif method is None or method == 'None':
        return df

    else:
        print(f"FillNA method `{str(method)}` not known, returning original")
        return df


text:
#! /usr/bin/env python
# -*- coding:utf-8
# K nearest neighbor
# Sample Data : https://archive.ics.uci.edu/ml/datasets/Breast+Cancer+Wisconsin+%28Original%29
# Special thanks: Harisson@pythonprogramming.net

'''
    File name: knn.py
    Author: chimney37
    Date created: 8/11/2017
    Python Version: 3.62
'''
import sys
import numpy as np
from sklearn import cross_validation, neighbors
import pandas as pd

# Column attributes (id,clump_thickness,uniform_cell_size,
# uniform_cell_shape,marginal_adhesion,
# single_epi_cell_size,bare_nuclei,bland_chromation,
# normal_nucleoli,mitoses,class)
# For more info, see Attribute information:
# https://archive.ics.uci.edu/ml/machine-learning-databases/breast-cancer-wisconsin/breast-cancer-wisconsin.names
df = pd.read_csv('breast-cancer-wisconsin.data.txt')
df.replace('?', -sys.maxsize, inplace=True)
df.drop(['id'], 1, inplace=True)

# drop the column(s) not used for training
X = np.array(df.drop(['class'], 1))

# labels (for prediction)
Y = np.array(df['class'])

# create training and testing samples
X_train, X_test, y_train, y_test = cross_validation.train_test_split(X, Y, test_size=0.2)

# Train the classifier
clf = neighbors.KNeighborsClassifier()
clf.fit(X_train, y_train)

accuracy = clf.score(X_test, y_test)
print('accuracy=', accuracy)

# crete example data for prediction
example_measures = np.array([[4, 2, 1, 1, 1, 2, 3, 2, 1], [4, 2, 1, 1, 1, 2, 3, 2, 1]])
example_measures = example_measures.reshape(len(example_measures), -1)

prediction = clf.predict(example_measures)
print(prediction)


text:
from __future__ import print_function

import numpy as np

def IsPowerOfTwo(i):
    """Returns true if all entries of i are powers of two, False otherwise.

    """
    return (i & (i - 1)) == 0 and i != 0

def Log2ofPowerof2(shape):
    """ Returns powers of two exponent for each element of shape

    """
    res = np.array(shape)
    for i in range(res.size):
        n = shape[i]
        assert (IsPowerOfTwo(n)), "Invalid input"
        ix = 0
        while n > 1:
            n //= 2
            ix += 1
        res[i] = ix
    return res

def Freq(i, N):
    """Outputs the absolute integers frequencies [0,1,...,N/2,N/2-1,...,1]

     in numpy fft convention as integer i runs from 0 to N-1.
     Inputs can be numpy arrays e.g. i (i1,i2,i3) with N (N1,N2,N3)
                                  or i (i1,i2,...) with N
     Both inputs must be integers.
     All entries of N must be even.

    """
    assert (np.all(N % 2 == 0)), "This routine only for even numbers of points"
    return i - 2 * (i >= (N / 2)) * (i % (N / 2))

def rfft2_reals(shape):
    """Pure reals modes in 2d rfft array. (from real map shape, not rfft array)

    """
    N0, N1 = shape
    fx = [0]
    fy = [0]
    if N1 % 2 == 0: fx.append(0); fy.append(N1 // 2)
    if N0 % 2 == 0: fx.append(N0 // 2); fy.append(0)
    if N1 % 2 == 0 and N0 % 2 == 0: fx.append(N0 // 2); fy.append(N1 // 2)
    return np.array(fx), np.array(fy)

def upgrade_map(LD_map, HD_res):
    """
    Upgrade LD_map to a higher resolution map, using rfft and back.
    :param LD_map: Must have shape entries powers of two.
    :param HD_res:
    :return: Same map at higher resolution.
    """
    LD_res = Log2ofPowerof2(LD_map.shape)
    if np.all(LD_res == HD_res): return LD_map
    assert np.all(HD_res >= LD_res)
    HD_rshape = (2 ** HD_res[0], 2 ** (HD_res[1] - 1) + 1)
    HD_shape = (2 ** HD_res[0], 2 ** HD_res[1])

    rfft = np.fft.rfft2(LD_map)
    ret_rfft = np.zeros(HD_rshape, dtype=complex)
    ret_rfft[0:rfft.shape[0] / 2 + 1, 0:rfft.shape[1]] = rfft[0:rfft.shape[0] / 2 + 1, :]  # positive frequencies
    # negative frequencies :
    ret_rfft[HD_rshape[0] - rfft.shape[0] + rfft.shape[0] / 2:, 0:rfft.shape[1]] = rfft[rfft.shape[0] / 2:, :]
    fac_LDrfft2HDrfft = 2 ** (HD_res[0] - LD_res[0] + HD_res[1] - LD_res[1])
    return np.fft.irfft2(ret_rfft, HD_shape) * fac_LDrfft2HDrfft

def subsample(HD_map, LD_res):
    """Simple subsampling of map.

    """
    HD_res = Log2ofPowerof2(HD_map.shape)
    if np.all(LD_res == HD_res): return HD_map.copy()
    assert np.all(HD_res >= LD_res)
    return HD_map[0::2 ** (HD_res[0] - LD_res[0]), 0::2 ** (HD_res[1] - LD_res[1])]

def supersample(LD_map, HD_shape):
    """Simple hypersampling of map.

    """
    if LD_map.shape == HD_shape: return LD_map.copy()
    assert np.all(np.array(HD_shape) > np.array(LD_map.shape))
    assert np.all(np.array(HD_shape) % np.array(LD_map.shape) == 0.)
    HDmap = np.zeros(HD_shape)
    fac0, fac1 = (HD_shape[0] / LD_map.shape[0], HD_shape[1] / LD_map.shape[1])
    for i in range(fac0):
        for j in range(fac1):
            HDmap[i::fac0, j::fac1] = LD_map
    return HDmap

def degrade(HD_map, LD_shape):
    if np.all(HD_map.shape <= LD_shape): return HD_map.copy()
    fac0, fac1 = (HD_map.shape[0] / LD_shape[0], HD_map.shape[1] / LD_shape[1])
    assert fac0 * LD_shape[0] == HD_map.shape[0] and fac1 * LD_shape[1] == HD_map.shape[1], (
        (fac0, fac1), LD_shape, HD_map.shape)

    ret = np.zeros(LD_shape, dtype=HD_map.dtype)
    for _i in range(fac0):
        sl0 = slice(_i, HD_map.shape[0], fac0)
        for _j in range(fac1):
            sl1 = slice(_j, HD_map.shape[1], fac1)
            ret += HD_map[sl0, sl1]
    return ret * (1. / (fac0 * fac1))

def degrade_mask(mask, LD_shape):
    # FIXME :
    dmask = degrade(mask, LD_shape)
    return dmask  # * (dmask >= 1.)

def udgrade_rfft2(rfft2map, shape, norm=False):
    assert norm == False, 'not implemented'  # norm. factor for rfft normalization maps. Here just shuffling indices.
    if shape == (rfft2map.shape[0], 2 * (rfft2map.shape[1] - 1)):
        return rfft2map
    assert np.all([s % 2 == 0 for s in shape]), shape
    assert rfft2map.shape[0] % 2 == 0
    rshape = np.array((shape[0], shape[1] / 2 + 1))
    if np.all(np.array(rshape) >= rfft2map.shape):
        return _upgrade_rfft2(rfft2map, shape)
    elif np.all(np.array(rshape) <= rfft2map.shape):
        return _degrade_rfft2(rfft2map, shape)
    else:
        assert 0, 'not implemented'

def _degrade_rfft2(rfft2map, LDshape):
    ret = np.zeros((LDshape[0], LDshape[0] // 2 + 1), dtype=complex)
    ret[0:LDshape[0] // 2 + 1, :] = rfft2map[0:LDshape[0] // 2 + 1, 0:ret.shape[1]]
    ret[LDshape[0] // 2::] = rfft2map[rfft2map.shape[0] - LDshape[0] // 2:, 0:ret.shape[1]]
    # Corrections for pure reals and (-k) = k* :
    ret[LDshape[0] // 2 + 1:, -1] = ret[1:LDshape[0] // 2, -1][::-1].conj()
    ret[rfft2_reals(LDshape)].imag = 0.
    return ret

def _upgrade_rfft2(rfft2map, HDshape):
    ret = np.zeros((HDshape[0], HDshape[0] / 2 + 1), dtype=complex)
    # positive 0axis frequencies : (including N/2 + 1, which is pure real.
    ret[0:rfft2map.shape[0] // 2 + 1, 0:rfft2map.shape[1]] = rfft2map[0:rfft2map.shape[0] // 2 + 1, 0:rfft2map.shape[1]]
    # Negative 0axis freq.
    ret[HDshape[0] - rfft2map.shape[0] // 2:HDshape[0], 0:rfft2map.shape[1]] \
        = rfft2map[rfft2map.shape[0] // 2:, 0:rfft2map.shape[1]]
    return ret

max_stars_repo_name: Sebastian-Belkner/LensIt

max_issues_repo_name: Sebastian-Belkner/LensIt

max_forks_repo_name: Sebastian-Belkner/LensIt


text:
import matplotlib
matplotlib.use('Agg')
import matplotlib.pyplot as plt
import gzip
import os
import pickle
from itertools import groupby, repeat
from subprocess import check_output, call
from os.path import dirname
import forgi.graph.bulge_graph as fgb
import numpy as np
from shutil import which, move
from collections import Counter
from os import remove
from tempfile import gettempdir
from math import ceil

def save_model(model, file_path):
    """ Save a pysster.Model object.

    This function creates two files: a pickled version of the pysster.Model object and
    an hdf5 file of the actual keras model (e.g. if file_path is 'model' two files are
    created: 'model' and 'model.h5')

    Parameters
    ----------
    model : pysster.Model
        A Model object.

    file_path : str
        A file name.
    """
    with gzip.open(file_path, "wb") as handle:
        pickle.dump(model.params, handle, pickle.HIGHEST_PROTOCOL)
    model.model.save("{}.h5".format(file_path))

def load_model(file_path, bluewhale=False):
    """ Load a pysster.Model object.

    Parameters
    ----------
    file_path : str
       A file containing a pickled pysster.Model object (file_path.h5 must also exist, see save_model()).

    Returns
    -------
    model : pysster.Model
        A Model object.
    """
    from pysster.Model import Model
    from keras.models import load_model as load_keras
    from janggu_layers import Reverse,Complement
    if not os.path.exists(file_path):
        raise RuntimeError("Path not found.")
    if not os.path.exists("{}.h5".format(file_path)):
        raise RuntimeError("HDF5 file not found.")
    with gzip.open(file_path, "rb") as handle:
        params = pickle.load(handle)
    custom_objects={'Reverse':Reverse,'Complement':Complement}
    model = Model(params, None, bluewhale=bluewhale)
    model.model = load_keras("{}.h5".format(file_path),custom_objects=custom_objects)
    return model

def save_data(data, file_path):
    """ Save a pysster.Data object.

    The object will be pickled to disk.

    Parameters
    ----------
    file_path : str
        A file name.
    """
    with gzip.open(file_path, "wb") as handle:
        pickle.dump(data, handle, pickle.HIGHEST_PROTOCOL)

def load_data(file_path):
    """ Load a pysster.Data object.

    Parameters
    ----------
    file_path : str
        A file containing a pickled pysster.Data object.

    Returns
    -------
    data : pysster.Data
        The Data object loaded from file.
    """
    with gzip.open(file_path, "rb") as handle:
        return pickle.load(handle)

def get_handle(file_name, mode):
    if file_name[-2:] == "gz":
        return gzip.open(file_name, mode)
    return open(file_name, mode)

def parse_fasta(handle, joiner = ""):
    delimiter = lambda line: line.startswith('>')
    for is_header, block in groupby(handle, delimiter):
        if is_header:
            header = next(block)[1:].rstrip()
        else:
            yield(header, joiner.join(line.rstrip() for line in block))

def annotate_structures(input_file, output_file):
    """ Annotate secondary structure predictions with structural contexts.

    Given dot-bracket strings this function will annote every character as
    either 'H' (hairpin), 'S' (stem), 'I' (internal loop/bulge), 'M' (multi loop), 'F' (5-prime)
    or 'T' (3-prime). The input file must be a fasta formatted file and each sequence and structure
    must span a single line:

    '>header
    'CCCCAUAGGGG
    '((((...)))) (-3.3)

    This is the default format of e.g. RNAfold. The output file will contain the annotated string:

    '>header
    'CCCCAUAGGGG
    'SSSSHHHSSSS

    Parameters
    ----------
    input_file : str
        A fasta file containing secondary structure predictions.

    output_file : str
        A fasta file with secondary structure annotations.
    """
    handle_in = get_handle(input_file, "rt")
    handle_out = get_handle(output_file, "wt")
    for header, entry in parse_fasta(handle_in, "_"):
        entry = entry.split("_")
        bg = fgb.BulgeGraph.from_dotbracket(entry[1].split()[0])
        handle_out.write(">{}\n".format(header))
        handle_out.write("{}\n{}\n".format(entry[0], bg.to_element_string().upper()))
    handle_in.close()
    handle_out.close()

def predict_structures(input_file, output_file, num_processes=None, annotate=False):
    """ Predict secondary structures for RNA sequences.

    This is a convenience function to get quick RNA secondary structure predictions. The function
    will try to use the ViennaRNA python bindings or the RNAfold binary to perform predictions. If
    neither can be found the function raises a RuntimeError. Using the ViennaRNA python bindings
    is preferred as it is faster.

    Entries of the output file look as follows if annotate = False:

    '>header
    'CCCCAUAGGGG
    '((((...)))) (-3.3)

    If annotate = True the annotated structure string instead of the dot-bracket string will be printed:

    '>header
    'CCCCAUAGGGG
    'SSSSHHHSSSS

    Have a look at the annotate_structures() function for more information about annotated structure strings.

    Warning: Due to the way Python works spinning up additional processes means copying the complete
    memory of the original process, i.e. if the original processes already uses 5 GB of RAM each additional
    process will use 5 GB as well.

    Parameters
    ----------
    input_file : str
        A fasta file with RNA sequences.

    output_file : str
        A fasta file with sequences and structures.

    num_processes : int
        The number of parallel processes to use for prediction. (default: number of available cores)

    annotate : bool
        Output the annotated structure string instead of the dot-bracket string. (default: false)
    """
    from multiprocessing import Pool

    try:
        from RNA import fold
        predictor = _predict_rnalib
    except:
        if which("RNAfold") == None:
            raise RuntimeError("Error: Neither ViennaRNA python bindings nor RNAfold executable found.")
        predictor = _predict_binary
    handle = get_handle(input_file, "rt")
    if num_processes == None:
        num_processes = max(1, int(os.cpu_count()/2))
    with Pool(num_processes) as pool:
        if annotate:
            data = pool.starmap(func = _predict_and_annotate,
                                iterable = zip(parse_fasta(handle), repeat(predictor)),
                                chunksize = 2)
            formatter = ">{}\n{}\n{}\n"
        else:
            data = list(pool.imap(func = predictor,
                                  iterable = parse_fasta(handle),
                                  chunksize = 2))
            formatter = ">{}\n{}\n{} ({})\n"
    handle.close()
    handle = get_handle(output_file, "wt")
    for entry in data:
        handle.write(formatter.format(*entry))
    handle.close()

def _predict_rnalib(fasta_entry):
    from RNA import fold
    return (*fasta_entry, *fold(fasta_entry[1]))

def _predict_binary(fasta_entry):
    out = check_output("echo {} | RNAfold --noPS".format(fasta_entry[1]), shell = True).split()
    return (*fasta_entry, out[1].decode("utf-8"), float(out[-1][:-1].replace(b"(", b"")))

def _predict_and_annotate(fasta_entry, predict_function):
    predict_entry = predict_function(fasta_entry)
    bg = fgb.BulgeGraph.from_dotbracket(predict_entry[2])
    return (predict_entry[0], predict_entry[1], bg.to_element_string().upper())

def auROC(labels, predictions):
    from sklearn.metrics import auc, roc_curve
    fpr, tpr, _ = roc_curve(labels, predictions)
    return fpr, tpr, auc(fpr, tpr)

def auPR(labels, predictions):
    from sklearn.metrics import precision_recall_curve, average_precision_score
    precision, recall, _ = precision_recall_curve(labels, predictions)
    return precision, recall, average_precision_score(labels, predictions)

def performance_report(labels, predictions):
    from sklearn.preprocessing import label_binarize
    from sklearn.metrics import precision_recall_fscore_support
    classes =  list(range(labels.shape[1]))
    roc_aucs, pr_aucs  = [], []
    if len(classes) == 2:
        roc_aucs = [auROC(labels[:, 0], predictions[:, 0])[2]] * 2
        pr_aucs = [auPR(labels[:, 0], predictions[:, 0])[2]] * 2
        labels = label_binarize(np.argmax(labels, axis = 1), classes = classes)
    else:
        for x in classes:
            roc_aucs.append(auROC(labels[:, x], predictions[:, x])[2])
            pr_aucs.append(auPR(labels[:, x], predictions[:, x])[2])
    if not np.isclose(np.sum(predictions, axis=1), 1).all():
        # multi-label classification
        y_pred = predictions > 0.5
        y_pred.dtype = np.uint8
    else:
        y_pred = label_binarize(np.argmax(predictions, axis = 1), classes = classes)
    prec_recall_f1_support = precision_recall_fscore_support(labels, y_pred)
    report = np.empty((len(classes), 6))
    for x in classes:
        report[x,:] = [prec_recall_f1_support[0][x], prec_recall_f1_support[1][x],
                       prec_recall_f1_support[2][x], roc_aucs[x],
                       pr_aucs[x], prec_recall_f1_support[3][x]]
    return report

def get_performance_report(labels, predictions):
    """ Get a performance overview of a classifier.

    The report contains precision, recall, f1-score, ROC-AUC and Precision-Recall-AUC for every
    class (in a 1 vs. all approach) and weighted averages (weighted by the the number
    of sequences 'n' in each class).

    Parameters
    ----------
    labels : numpy.ndarray
        A binary matrix of shape (num sequences, num classes) containing the true labels.

    predictions : numpy.ndarray
        A matrix of shape (num sequences, num classes) containing predicted probabilites.

    Returns
    -------
    report : str
        Summary table of the above mentioned performance measurements.
    """
    classes =  list(range(labels.shape[1]))
    report = performance_report(labels, predictions)
    out = []
    out.append("             precision    recall  f1-score   roc-auc    pr-auc          n")
    formatter = "{:>12}" + "{:>10.3f}" * 5 + "  |" + "{:>8}"
    for x in classes:
        out.append(formatter.format("class_{}".format(x), *report[x, 0:-1], int(report[x,-1])))
    out.append('\n')
    out.append(formatter.format(
        "weighted avg",
        *np.sum(report[:,0:-1] * report[:,-1, np.newaxis], axis=0)/int(sum(report[:,-1])),
        " "
    ))
    out.append('\n')
    return '\n'.join(out)

def plot_roc(labels, predictions, file_path):
    """ Get ROC curves for every class.

    In the case of more than two classes the comparisons will be performed in a 1 vs. all
    approach (i.e. you get one curve per class).

    Parameters
    ----------
    labels : numpy.ndarray
        A binary matrix of shape (num sequences, num classes) containing the true labels.

    predictions : numpy.ndarray
        A matrix of shape (num sequences, num classes) containing predicted probabilites.

    file_path : str
        The file the plot should be saved to.
    """
    classes = list(range(labels.shape[1]))
    colors = _get_colors(len(classes))
    fig, ax = plt.subplots(nrows = 1, ncols = 1, figsize = (4.6666,4))
    _hide_top_right(ax)
    ax.plot([0, 1], [0, 1], color = 'black', linewidth = 1, linestyle = '--')
    ax.set_xlabel('False Positive Rate')
    ax.set_ylabel('True Positive Rate')
    if len(classes) == 2:
        fpr, tpr, roc_auc = auROC(labels[:, 0], predictions[:, 0])
        label = 'AUC = {:.3f}'.format(roc_auc)
        ax.plot(fpr, tpr, linewidth = 2.2, color = colors[0], label = label)
    else:
        for x in classes:
            fpr, tpr, roc_auc = auROC(labels[:, x], predictions[:, x])
            label = 'AUC class_{} = {:.3f}'.format(x, roc_auc)
            ax.plot(fpr, tpr, linewidth = 2.2, color = colors[x], label = label)
    ax.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0., framealpha=1)
    fig.savefig(file_path, bbox_inches = 'tight')
    plt.close(fig)

def plot_prec_recall(labels, predictions, file_path):
    """ Get Precision-Recall curves for every class.

    In the case of more than two classes the comparisons will be performed in a 1 vs. rest
    approach (i.e. you get one curve per class).

    Parameters
    ----------
    labels : numpy.ndarray
        A binary matrix of shape (num sequences, num classes) containing the true labels.

    predictions : numpy.ndarray
        A matrix of shape (num sequences, num classes) containing predicted probabilites.

    file_path : str
        The file the plot should be saved to.
    """
    classes = list(range(labels.shape[1]))
    colors = _get_colors(len(classes))
    fig, ax = plt.subplots(nrows = 1, ncols = 1, figsize = (4.6666,4))
    _hide_top_right(ax)
    ax.plot([0, 1], [0, 1], color = 'white', linewidth = 1, linestyle = '--')
    ax.set_xlabel('Recall')
    ax.set_ylabel('Precision')
    if len(classes) == 2:
        precision, recall, prec_auc = auPR(labels[:,0], predictions[:, 0])
        label = 'AUC = {:.3f}'.format(prec_auc)
        ax.plot(recall, precision, linewidth = 2.2, color = colors[0], label = label)
    else:
        for x in classes:
            precision, recall, prec_auc = auPR(labels[:, x], predictions[:, x])
            label = 'AUC class_{} = {:.3f}'.format(x, prec_auc)
            ax.plot(recall, precision, linewidth = 2.2, color = colors[x], label = label)
    ax.legend(bbox_to_anchor=(1.05, 1), loc=2, borderaxespad=0., framealpha=1)
    fig.savefig(file_path, bbox_inches = 'tight')
    plt.close(fig)

def save_as_meme(logos, file_path):
    """ Save sequence (or structure) motifs in MEME format.

    Parameters
    ----------
    logos : [pysster.Motif]
        A list of Motif objects.

    file_path : str
        The name of the output text file.
    """
    alphabet = logos[0].alphabet
    with open(file_path, "wt") as handle:
        handle.write("MEME version 4\n\nALPHABET= {}\n\nstrands: + -\n\n".format(alphabet))
        handle.write("Background letter frequencies (from uniform background):\n")
        header = ""
        for c in alphabet:
            header += "{} {:7.5f} ".format(c, 1/len(alphabet))
        handle.write(header[:-1] + '\n')
        for i, logo in enumerate(logos):
            pwm = logo.pwm
            handle.write("\nMOTIF motif_{} motif_{}\n\n".format(i,i))
            handle.write("letter-probability matrix: alength= {} w= {} nsites= 20 E= 0\n".format(
                len(alphabet),
                pwm.shape[0]
            ))
            for row in range(pwm.shape[0]):
                handle.write("  {}\t\n".format("\t  ".join(str(round(x, 6)) for x in pwm[row,:])))

def run_tomtom(motif_file, output_folder, database, options = None):
    """ Compare a MEME file against a database using TomTom.

    Default options string: "-min-overlap 5 -verbosity 1 -xalph -evalue -thresh 0.1"

    Parameters
    ----------
    motif_file : str
        A MEME file.

    output_folder : str
        The folder the TomTom output will be saved in.

    database : str
        A MEME file serving as the database to compare against.

    option : str
        Command line options passed to TomTom.
    """
    if which("tomtom") == None:
        raise RuntimeError("Error: tomtom executable not found.")
    if output_folder[-1] != "/":
        output_folder += "/"
    if not os.path.isdir(output_folder):
        os.makedirs(output_folder)
    if options == None:
        options = "-min-overlap 5 -verbosity 1 -xalph -evalue -thresh 0.1"
    command = "tomtom {}".format(options)
    call("{} -oc {} {} {}".format(command, output_folder, motif_file, database), shell = True)

def softmax(x):
    x = np.exp(x - np.max(x))
    return x / x.sum(axis = 0)

def _hide_top_right(ax):
    ax.spines['right'].set_visible(False)
    ax.spines['top'].set_visible(False)
    ax.yaxis.set_ticks_position('left')
    ax.xaxis.set_ticks_position('bottom')

def plot_motif_summary(position_max, mean_acts, kernel, file_path):
    from PIL import Image
    classes = []
    ylim_hist, ylim_mean = 0, 0
    for i, hist in enumerate(position_max):
        if len(hist) == 0:
            print("Warning: class {} did not activate kernel {}. No plots were created.".format(
                i, kernel
            ))
        else:
            classes.append(i)
            ylim_hist = max(ylim_hist, Counter(hist).most_common(1)[0][1])
            ylim_mean = max(ylim_mean, max(mean_acts[i][0] + mean_acts[i][1]))
    xlim = len(mean_acts[classes[0]][0]) + 1
    old_fontsize = matplotlib.rcParams['font.size']
    matplotlib.rcParams.update({'font.size': 30})
    files = []
    n_per_plot = 3
    n_plots = ceil(len(classes)/n_per_plot)
    class_idx = -1
    for plot_id in range(n_plots):
        classes_left = len(classes) - plot_id*n_per_plot
        classes_this_plot = min(n_per_plot, classes_left)
        fig, ax = plt.subplots(nrows = 2,
                               ncols = classes_this_plot,
                               figsize = (19*classes_this_plot, 12))
        for class_num in range(classes_this_plot):
            class_idx += 1
            # histograms
            ax.flat[class_num].set_ylim((0, ylim_hist))
            ax.flat[class_num].hist(position_max[classes[class_idx]], histtype="stepfilled",
                                    bins = xlim, range = (0, xlim))
            ax.flat[class_num].set_xlabel("sequence position")
            ax.flat[class_num].set_ylabel("counts")
            ax.flat[class_num].set_title("kernel {}, class_{}, (n = {})".format(
            kernel, classes[class_idx], len(position_max[classes[class_idx]])
            ))
            _hide_top_right(ax.flat[class_num])
            # mean activations
            ax.flat[class_num + classes_this_plot].set_ylim((0, ylim_mean))
            ax.flat[class_num + classes_this_plot].fill_between(list(range(1, xlim)),
                                                                mean_acts[classes[class_idx]][0] - mean_acts[classes[class_idx]][1],
                                                                mean_acts[classes[class_idx]][0] + mean_acts[classes[class_idx]][1],
                                                                alpha = 0.1)
            ax.flat[class_num + classes_this_plot].plot(list(range(1, xlim)),
                                                        mean_acts[classes[class_idx]][0],
                                                        linewidth = 5.0)
            ax.flat[class_num + classes_this_plot].set_xlabel("sequence position")
            ax.flat[class_num + classes_this_plot].set_ylabel("activation")
            _hide_top_right(ax.flat[class_num + classes_this_plot])
        fig.tight_layout()
        files.append("{}/plotsum{}.png".format(gettempdir(), plot_id))
        fig.savefig(files[-1])
        plt.close(fig) # fig.clf() before close() seems to release memory faster
    matplotlib.rcParams.update({'font.size': old_fontsize})
    if len(files) == 1:
        move(files[0], file_path)
    else:
        images = []
        for file_name in files:
            images.append(Image.open(file_name))
        combine_images(images, file_path)
        for file_name in files:
            remove(file_name)

def plot_violins(data, kernel, file_path):
    old_fontsize = matplotlib.rcParams['font.size']
    matplotlib.rcParams.update({'font.size': 15})
    num_plots = len(data)
    labels = ["class_{}".format(x) for x in range(num_plots)]
    fig, ax = plt.subplots(nrows = 1, ncols = 1, figsize = (max(5, num_plots), 5))
    _hide_top_right(ax)
    ax.grid(axis = "y", alpha = 0.3)
    ax.set_title("Activations, kernel {}".format(kernel))
    ax.set_ylabel("max activations")
    parts = ax.violinplot(data, showmeans = True, showextrema = True)
    ax.set_ylim(bottom = 0)
    parts['cmeans']._linewidths = [2]
    parts['cmins']._linewidths = [2]
    parts['cmaxes']._linewidths = [2]
    parts['cbars']._linewidths = [2]
    ax.get_xaxis().set_tick_params(length=0)
    ax.set_xticks(np.arange(1, len(labels) + 1))
    ax.set_xticklabels(labels)
    ax.set_xlim(0.25, len(labels) + 0.75)
    plt.xticks(rotation=90)
    fig.subplots_adjust(bottom=0.2)
    plt.tight_layout()
    fig.savefig(file_path)
    plt.close(fig)
    matplotlib.rcParams.update({'font.size': old_fontsize})

def plot_motif(logo, file_path, colors_sequence, colors_structure):
    from PIL import Image
    if isinstance(logo, tuple):
        img1 = logo[0].plot(colors_sequence, scale=0.75)
        img2 = logo[1].plot(colors_structure, scale=0.75)
        img = Image.new("RGB", (img1.size[0], img1.size[1]+img2.size[1]))
        img.paste(img1, (0, 0))
        img.paste(img2, (0, img1.size[1]))
        img1.close()
        img2.close()
    else:
        img = logo.plot(colors_sequence, scale=0.75)
    img.save(file_path)
    img.close()
    plt.close('all')

def _set_sns_context(n_kernel):
    import seaborn as sns
    if n_kernel <= 25:
        sns.set_context("notebook", rc={"ytick.labelsize":26})
    elif 25 < n_kernel <= 50:
        sns.set_context("notebook", rc={"ytick.labelsize":22})
    elif 50 < n_kernel <= 75:
        sns.set_context("notebook", rc={"ytick.labelsize":14})
    elif 75 < n_kernel <= 100:
        sns.set_context("notebook", rc={"ytick.labelsize":8})
    else:
        sns.set_context("notebook", rc={"ytick.labelsize":5})

def _get_colors(x):
    import seaborn as sns
    palette = ["hls", "Set1"][x < 10]
    return sns.color_palette(palette, x, 0.6)

def _plot_heatmap(file_path, data, class_id, classes = None):
    import seaborn as sns
    _set_sns_context(data.shape[1])
    n_classes = len(set(class_id))
    palette = _get_colors(n_classes)
    colors = [palette[x] for x in class_id]
    g = sns.clustermap(data = data.T, method = "ward", metric = "euclidean",
                       cmap = "RdBu_r", xticklabels = False, yticklabels = True,
                       figsize = (30,25), row_cluster = True, col_cluster = True,
                       linewidths = 0, col_colors = colors, robust = True,
                       z_score = 0, cbar_kws={"ticks":[-1.5,0,+1.5]})
    g.ax_col_dendrogram.set_xlim([0,1e-10])
    g.ax_col_dendrogram.set_ylim([0,1e-10])
    plt.setp(g.ax_heatmap.get_yticklabels(), rotation=0)
    sns.set(font_scale=2.8)
    if classes == None:
        classes = list(range(n_classes))
    for x in range(n_classes):
        g.ax_col_dendrogram.bar(0, 0, color=palette[x],
                                label="class_{}".format(classes[x]), linewidth=0)
    g.ax_col_dendrogram.legend(loc = "center", ncol = min(6, n_classes))
    cax = plt.gcf().axes[-1]
    cax.tick_params(labelsize=25)
    plt.savefig(file_path, bbox_inches = 'tight')
    plt.close('all')
    sns.reset_orig()

def combine_images(images, output_file):
    from PIL import Image
    widths, heights = zip(*(i.size for i in images))
    new_im = Image.new('RGB', (max(widths), sum(heights)), "#ffffff")
    y_offset = 0
    for im in images:
        new_im.paste(im, (0, y_offset))
        im.close()
        y_offset += im.size[1]
    new_im.save(output_file)
    new_im.close()
    plt.close('all')

# helper function; np.argmax always returns the first occurrence of the max value
# even if it occurs multiple times; this function randomly selects one of those instead
def randargmax(data):
    rtol, atol = 1e-09, 0.0
    result = np.empty((data.shape[0],), dtype=np.uint32)
    max_val = np.max(data, axis=1)
    for x in range(data.shape[0]):
        result[x] = np.random.choice(
            np.where(abs(data[x, ] - max_val[x]) <= np.maximum(rtol * np.maximum(abs(data[x, ]), abs(max_val[x])), atol))[0])
    return result

def html_report(sorted_idx, scores, folder, class_num, size=None):
    handle = open(folder+"summary.html", "wt")
    handle.write('<!doctype html>\n<html>\n<head>\n<meta charset="UTF-8">\n')
    handle.write('<title>Kernel Summary</title>\n<style media="screen" type="text/css">\n')
    handle.write('#report {white-space: nowrap;}\n')
    handle.write('td {text-align: center; font-weight: bold; padding: 20px;}\n')
    handle.write('table {margin: 0 auto; border-collapse: collapse;}\n')
    handle.write('tr:nth-child(even) {background-color: #f2f9ff;}</style>\n</head>\n')
    handle.write('<body>\n<div id="report">\n<table>\n')
    for kernel in sorted_idx:
        handle.write('<tr>\n<td>Kernel {}<br/>score = {:.3f}</td>\n'.format(kernel, scores[kernel]))
        if size == None:
            handle.write('<td><img src="motif_kernel_{}.png" height=150/></td>\n'.format(kernel))
        else:
            handle.write('<td><img src="motif_kernel_{}.png" height=150/><br>'.format(kernel))
            handle.write('<img src="additional_features_kernel_{}.png" height={}/></td>\n'.format(kernel, size))
        handle.write('<td><img src="activations_kernel_{}.png" height=300/></td>\n'.format(kernel))
        handle.write('<td><img src="position_kernel_{}.png" height={}/></td>\n</tr>\n'.format(
            kernel, 300 * max(ceil(class_num/3), 1)))
    handle.write('</table>\n</div>\n</body>\n</html>')
    handle.close()

def plot_positionwise(add_data, identifiers, file_path):
    old_fontsize = matplotlib.rcParams['font.size']
    matplotlib.rcParams.update({"font.size": 30})
    x = list(range(1, add_data[0][0].shape[0]+1))
    fig, ax = plt.subplots(nrows=len(add_data), ncols=1,
                           figsize=(int(26*(len(x)/20)), 7*len(add_data)))
    if not isinstance(ax, np.ndarray):
        ax = [ax]
    for i in range(len(add_data)):
        mean, std = np.mean(add_data[i], axis=0), np.std(add_data[i], axis=0)
        ax[i].fill_between(x, mean-std, mean+std, color='orange', alpha = 0.1)
        ax[i].plot(x, mean, 'o-', linewidth=5.0, markersize=15.0, color="orange")
        if len(x) <= 30:
            ax[i].set_xticks(x)
        elif len(x) > 75:
            ax[i].set_xticks(list(range(1, max(x), 3)))
        else:
            ax[i].set_xticks(list(range(1, max(x), 2)))
        ax[i].set_ylabel("{}".format(identifiers[i]))
        _hide_top_right(ax[i])
    plt.tight_layout()
    fig.savefig(file_path)
    plt.close(fig)
    matplotlib.rcParams.update({"font.size": old_fontsize})

max_stars_repo_name: sproft/pysster

max_issues_repo_name: sproft/pysster

max_forks_repo_name: sproft/pysster


text:
from c0106_find_records import find_records
from c0112_define_records import define_records
from c0113_format_source import format_source
from c0114_format_truncate import format_truncate
from c0116_coregister_formatted import coregister_formatted
from c0117_clean_save import clean_save
from c0118_segment_formatted import segment_formatted

import os
import matplotlib
import matplotlib.pyplot as plt
import numpy as np
import pandas as pd

def clean_data():
    """
    Clean the data so each record only includes information measured from the participant.

    If the wearable is taken off and then turned off - the extra time off-wrist is truncated
    If the wearable is left on between participants, so two records are combined in one - the record is broken into two
    If two wearables were used to (one on the left and one on the right wrist) - the data is coregistered.
    Timestamps are assigned to the measurements.
    Records are broken into segment for comparision - preStudy, Study, postStudy

    Here are some definitions:

    1. source - a description of data, uneditted measurements downloaded from the cloud client
    2. truncate - a description of data, source data only when on the participant's wrist

    3. timestamp - action taken on the data, create the lit of unix and time in minutes corresponding to the measurements
    4. coregister - action taken on the data, pair up records, for some studies two wearables were used for the same patient
    5. segment - action taken on the data, break the entire truncated record into a preStudy, Study, and postStudy period

    6. define record - processing step, establish the beginning and end of the record in unix time
    7. embedded record - type of record, a record that contains two separate recorded session from what could be two people
    embeddd records need to be broken up into two separate records, which is handled in this cleaning process

    How is the end of the record found?
    Use the temperature sensor. The room is significantly colder than the participant's wrist.
    If the temperature drops more than 2 deg C in 3 seconds - the wearable is assummed removed.
    If the temperature rises more than 2 deg C in 3 seconds - the wearable is assummed put on - only used for embedded records.
    """

    print("begin cleaning data")

    # find and list all records
    find_records()

    # define record - find record begin, end, length
    define_records()

    # timestamp the source
    format_source()

    # timestamp the source
    format_truncate()

    # coregister wearables
    coregister_formatted()

    # save processed data as clean
    clean_save()

    # segment clean records
    format_type = 'clean'
    segment_formatted(format_type)

    print("completed cleaning data")

max_stars_repo_name: jesnyder/MeasuredStress

max_issues_repo_name: jesnyder/MeasuredStress

max_forks_repo_name: jesnyder/MeasuredStress


text:
import numpy as np
import matplotlib as plt
import fastai as fa

x = np.random.rand(3,4 ,2,2)
print(x)

y = np.ones([3, 4], dtype=float)

print(y)

max_stars_repo_name: waeschinger/Deepsy

max_issues_repo_name: waeschinger/Deepsy

max_forks_repo_name: waeschinger/Deepsy


text:
from numbers import Number
import numpy as np

class ViewBase:
    def __init__(self, codes, lengths=None):
        if lengths is None:
            self._codes = codes.view(np.int32)
        else:
            starts = np.asanyarray(codes, dtype=np.int32)
            lengths = np.asanyarray(lengths, dtype=np.int32)
            if not lengths.size:
                self._codes = np.array([], dtype=np.int32)
            else:
                self._codes = np.hstack((starts[:, None], lengths[:, None])).flatten()

    def __eq__(self, other):
       return np.all(self._codes==other._codes)

    def __repr__(self):
        return f"{self.__class__.__name__}({self.starts}, {self.lengths})"

    @property
    def lengths(self):
        """The row lengths"""
        return self._codes[1::2]

    @property
    def starts(self):
        """The start index of each row"""
        return self._codes[::2]

    @property
    def ends(self):
        """The end index of each row"""
        return self.starts+self.lengths

    @property
    def n_rows(self):
        """Number of rows"""
        if isinstance(self.starts, Number):
            return 1
        return self.starts.size

    def empty_rows_removed(self):
        """Check wheter the `View` with certainty have no empty rows

        Returns
        -------
        bool
            Whether or not it is cerain that this view contins no empty rows
        """
        return hasattr(self, "empty_removed") and self.empty_removed

    def ravel_multi_index(self, indices):
        """Return the flattened indices of a set of array indices

        Parameters
        ----------
        indices : tuple
            Tuple containing the row- and column indices to ravel

        Returns
        -------
        array
            array containing the flattenened indices
        """
        return self.starts[indices[0]]+np.asanyarray(indices[1], dtype=np.int32)

    def unravel_multi_index(self, flat_indices):
        """Return array indices for a set of flat indices

        Parameters
        ----------
        indices : index_like
            flat indices to unravel

        Returns
        -------
        tuple
            tuple containing the unravelled row- and column indices
        """
        starts = self.starts
        rows = np.searchsorted(starts, flat_indices, side="right")-1
        cols = flat_indices-starts[rows]
        return rows, cols

    def index_array(self):
        """Return an array of broadcasted row indices"""
        diffs = np.zeros(self.size, dtype=np.int32)
        diffs[self.starts[1:]] = 1
        return np.cumsum(diffs)

class RaggedRow:
    def __init__(self, code):
        code =  np.atleast_1d(code).view(np.int32)
        self.starts = code[0]
        self.legths = code[1]
        self.ends = code[0]+code[1]

class RaggedShape(ViewBase):
    """ Class that represents the shape of a ragged array.

    Represents the same information as a list of row lengths.

    Parameters
    ----------
    codes : list or array_like
        Either a list of row lengths, or if ``is_coded=True`` an  array containing row-starts
        and row-lengths as 32-bit numbers.
    is_coded : bool, default=False
        if `False`, the `codes` are interpreted as row lengths.

    Attributes
    ----------
    starts
    lengths
    ends
    """
    def __init__(self, codes, is_coded=False):
        if is_coded: # isinstance(codes, np.ndarray) and codes.dtype==np.uint64:
            super().__init__(codes)
            self._is_coded = True
        else:
            lengths = np.asanyarray(codes, dtype=np.int32)
            starts = np.insert(lengths.cumsum(dtype=np.int32)[:-1], 0, np.int32(0))
            super().__init__(starts, lengths)
            self._is_coded = True

    def __repr__(self):
        return f"{self.__class__.__name__}({self.lengths})"

    def __str__(self):
        return str(self.lengths)

    def __getitem__(self, index):
        if not isinstance(index, slice) or isinstance(index, Number):
            return NotImplemented
        if isinstance(index, Number):
            index = [index]
        new_codes = self._codes.view(np.uint64)[index].copy().view(np.int32)
        new_codes[::2] -= new_codes[0]
        return self.__class__(new_codes, is_coded=True)

    @property
    def size(self):
        """The sum of the row lengths"""
        if not self.n_rows:
            return 0
        return self.starts[-1]+self.lengths[-1]

    def view(self, indices):
        """Return a view of a subset of rows

        Return a view with row information for the row given by `indices`

        Parameters
        ----------
        indices : index_like
            Used to index the rows

        Returns
        -------
        RaggedView
            RaggedView containing information to find the rows specified by `indices`
        """
        if isinstance(indices, Number):
            return RaggedRow(self._codes.view(np.uint64)[indices])
        return RaggedView(self._codes.view(np.uint64)[indices])

    def to_dict(self):
        """Return a `dict` of all necessary variables"""
        return {"codes": self._codes}

    @classmethod
    def from_dict(cls, d):
        """Load a `Shape` object from a dict of necessary variables

        Paramters
        ---------
        d : dict
            `dict` containing all the variables needed to initialize a RaggedShape

        Returns
        -------
        RaggedShape
        """
        if "offsets" in d:
            return cls(np.diff(d["offsets"]))
        else:
            return cls(d["codes"])

    @classmethod
    def asshape(cls, shape):
        """Create a `Shape` from either a list of row lengths or a `Shape`

        If `shape` is already a `RaggedShape`, do nothing. Else construct a new
        `RaggedShape` object

        Parameters
        ----------
        shape : RaggedShape or array_like

        Returns
        -------
        RaggedShape
        """
        if isinstance(shape, RaggedShape):
            return shape
        return cls(shape)

    def broadcast_values(self, values, dtype=None):
        """Broadcast the values in a column vector to the data of a ragged array

        The resulting array is such that a `RaggedArray` with `self` as shape will
        have the rows filled with the values in `values. I.e.
        ``RaggedArray(ret, self)[row, j] = values[row, 1]``

        Parameters
        ----------
        values : array_like
            column vectors with values to be broadcasted

        Returns
        -------
        array
            flat array with broadcasted values
        """
        values = np.asanyarray(values)
        assert values.shape == (self.n_rows, 1), (values.shape, (self.n_rows, 1))
        if self.empty_rows_removed():
            return self._broadcast_values_fast(values, dtype)
        values = values.ravel()
        broadcast_builder = np.zeros(self.size+1, dtype=dtype)
        broadcast_builder[self.ends[::-1]] -= values[::-1]
        broadcast_builder[0] = 0
        broadcast_builder[self.starts] += values
        func = np.logical_xor if values.dtype==bool else np.add
        return func.accumulate(broadcast_builder[:-1])

    def _broadcast_values_fast(self, values, dtype=None):
        values = values.ravel()
        broadcast_builder = np.zeros(self.size, dtype=dtype)
        broadcast_builder[self.starts[1:]] = np.diff(values)
        broadcast_builder[0] = values[0]
        func = np.logical_xor if values.dtype==bool else np.add
        func.accumulate(broadcast_builder, out=broadcast_builder)
        return broadcast_builder

class RaggedView(ViewBase):
    """Class to represent a view onto subsets of rows

    Same as RaggedShape, except without the constraint that the rows
    fill the whole data array. I.e. ``np.all(self.ends[:-1]==self.starts[1:])``
    does not necessarilty hold.

    Parameters
    ----------
    codes : array_like
        Either a list of row starts, or if `lengths` is provided an  array containing row-starts
        and row-lengths as 32-bit numbers.
    lengths : array_like, optional
        the lengths of the rows

    Attributes
    ----------
    starts
    lengths
    ends
    """

    def __getitem__(self, index):
        if isinstance(index, Number):
            return RaggedRow(self._codes.view(np.uint64)[index])

        return self.__class__(self._codes.view(np.uint64)[index])

    def get_shape(self):
        """ Return the shape of a ragged array containing the view's rows

        Returns
        -------
        RaggedShape
            The shape of a ragged array consisting of the rows in this view
        """
        if not self.n_rows:
            return RaggedShape(self._codes, is_coded=True)

        codes = self._codes.copy()
        np.cumsum(codes[1:-1:2], out=codes[2::2])
        codes[0] = 0
        return RaggedShape(codes, is_coded=True)

    def get_flat_indices(self):
        """Return the indices into a flattened array

        Return the indices of all the elements in all the
        rows in this view

        Returns
        -------
        array
        """
        if not self.n_rows:
            return np.ones(0, dtype=np.int32), self.get_shape()

        if self.empty_rows_removed():
            return self._get_flat_indices_fast()
        shape = self.get_shape()
        index_builder = np.ones(shape.size+1, dtype=np.int32)
        index_builder[shape.ends[::-1]] = 1-self.ends[::-1]
        index_builder[0] = 0
        index_builder[shape.starts] += self.starts
        np.cumsum(index_builder, out=index_builder)
        return index_builder[:-1], shape

    def _get_flat_indices_fast(self):
        shape = self.get_shape()
        index_builder = np.ones(shape.size, dtype=np.int32)
        index_builder[shape.starts[1:]] = np.diff(self.starts)-self.lengths[:-1]+1
        index_builder[0] = self.starts[0]
        np.cumsum(index_builder, out=index_builder)
        shape.empty_removed = True
        return index_builder, shape

max_stars_repo_name: jorgenwh/npstructures

max_issues_repo_name: jorgenwh/npstructures

max_forks_repo_name: jorgenwh/npstructures


text:
from train import Agent, utils as U
import gym
import tensorflow as tf
import numpy as np
import matplotlib.pyplot as plt

K = tf.keras
L = K.layers

class PPO(Agent):

    def __init__(self,
                 policy,
                 critic,
                 epsilon=.1,
                 T=16,
                 epochs=4,
                 batch_size=None,
                 optimizer=None,
                 transitions=-1,
                 **kwargs):
        super(PPO, self).__init__(transitions=transitions, **kwargs)
        self.policy = policy
        self.old = K.models.clone_model(policy)
        self.old.trainable = False
        self.critic = critic
        self.epsilon = epsilon
        self.T = T
        self.epochs = epochs
        self.batch_size = batch_size or T // 4
        self.optimizer = optimizer or K.optimizers.Adam()

    def act(self, state):
        probs = self.policy.predict(state[None])[0]
        action = np.random.choice(len(probs), p=probs)
        return action

    def on_step_end(self):
        if not self.old.built:
            S = self.transitions.sample(1).state
            self.old.predict(S)  # initialize weights
            self.update_old()
        if len(self.transitions) == self.T:
            self._learn()
            self.update_old()

    def update_old(self):
        self.old.set_weights(self.policy.get_weights())

    def _learn(self):
        data = self.transitions.get()
        self.transitions.reset()
        self.learn(data)

    def learn(self, data):
        S, A, R, Snext, dones = data
        T = len(R)
        batch_shape = (T, )
        gamma, lambd = self.gamma, self.lambd
        policy, old, critic = self.policy, self.old, self.critic

        old_probs = old.predict(S)
        old_probs = old_probs[range(old_probs.shape[0]), A]
        U.check_shape(old_probs, batch_shape)

        targets, deltas = self.compute_td_zero(data, V=critic.predict)
        advantages = self.compute_gae(deltas=deltas, dones=dones)

        indices = np.arange(T)
        for _ in range(self.epochs):
            np.random.shuffle(indices)
            for batch in np.array_split(indices, self.batch_size):
                self.optimize(S=S[batch],
                              A=A[batch],
                              old_probs=old_probs[batch],
                              advantages=advantages[batch],
                              targets=targets[batch])

    def optimize(self, S, A, old_probs, advantages, targets):
        batch_size = len(A)
        batch_shape = (batch_size, )
        epsilon, policy, critic = self.epsilon, self.policy, self.critic
        with tf.GradientTape() as tape:
            # Policy Objective
            probs = policy(S)
            probs = tf.gather(probs, A.reshape([-1, 1]), batch_dims=1)
            probs = tf.reshape(probs, [-1])
            U.check_shape(probs, batch_shape)
            ratios = probs / old_probs
            Lcpi = ratios * advantages
            Lclip = tf.clip_by_value(ratios, 1 - epsilon,
                                     1 + epsilon) * advantages
            policy_objective = tf.minimum(Lcpi, Lclip)
            U.check_shape(policy_objective, batch_shape)
            policy_objective = tf.reduce_mean(policy_objective)
            U.check_shape(policy_objective, ())
            # Critic Loss
            V = critic(S)
            V = tf.reshape(V, [-1])
            U.check_shape(V, batch_shape)
            critic_loss = tf.reduce_mean(tf.square(targets - V))
            U.check_shape(critic_loss, ())
            # Total Loss
            loss = -policy_objective + critic_loss
        grads = tape.gradient(loss, self.parameters)
        self.optimizer.apply_gradients(zip(grads, self.parameters))

    @property
    def parameters(self):
        if self._parameters is None:
            params = self.policy.trainable_variables + self.critic.trainable_variables
            params = U.unique(params)
            self._parameters = params
        return self._parameters

def main():
    env = gym.make('CartPole-v0')
    n_actions = env.action_space.n
    hidden = 16
    policy = K.Sequential([
        L.Dense(hidden, activation='relu'),
        L.Dense(n_actions, activation='softmax'),
    ])
    critic = K.Sequential([
        L.Dense(hidden, activation='relu'),
        L.Dense(1),
    ])
    agent = PPO(policy=policy, critic=critic, env=env)
    scores = agent.train(episodes=200)
    plt.plot(scores)
    plt.show()

if __name__ == '__main__':
    main()

max_stars_repo_name: marella/train

max_issues_repo_name: marella/train

max_forks_repo_name: marella/train


text:
import cv2
import numpy as np
from keras.models import load_model
from keras.preprocessing.image import img_to_array

from smile_module.base import SmileModel

class LenetSmile(SmileModel):
    name = 'LenetSmile'

    def _load_model(self) -> None:
        self.model = load_model(self.model_path)

    @staticmethod
    def _preprocess(cropped_face: np.ndarray) -> np.ndarray:
        def extract_roi(img, size):
            roi = cv2.resize(img, (size, size))
            roi = roi.astype("float") / 255.0
            roi = img_to_array(roi)
            return np.expand_dims(roi, axis=0)

        gray = cv2.cvtColor(cropped_face, cv2.COLOR_BGR2GRAY)
        region_of_interest = extract_roi(gray, size=28)
        return region_of_interest

    def _get_output(self, inp: np.ndarray) -> float:
        (_, smile) = self.model.predict(inp)[0]
        return smile

    def compute_score(self, cropped_face: np.ndarray) -> float:
        inp = self._preprocess(cropped_face)
        return self._get_output(inp)