text:
import tensorflow as tf
import time
import numpy as np
from sklearn.metrics import precision_score, recall_score, f1_score

class TextClassifier(object):
    """模型超参数"""

    def __init__(self,
                 n_hidden=128,
                 num_classes=2,
                 epoch=20,
                 max_seq_len=300,
                 embedding_size=50,
                 learning_rate=0.0005,
                 vocab_size=None,
                 num_layers=1,
                 wordVector=None,
                 preEmbedding=True,
                 isAttention=False,
                 filter_sizes=None,
                 num_filters=256,
                 num_blocks=2,
                 head=8):
        if filter_sizes is None:
            filter_sizes = [1, 2, 3, 4]
        self.n_hidden = n_hidden  # LSTM单元隐藏层神经元个数
        self.num_classes = num_classes  # 标签个数
        self.epoch = epoch  # 训练次数
        self.max_seq_len = max_seq_len  # 序列最大长度
        self.embedding_size = embedding_size  # 词向量维度
        self.learning_rate = learning_rate  # 学习率
        self.vocab_size = vocab_size  # 词典大小，embedding嵌入用
        self.num_layers = num_layers  # LSTM网络层数
        self.isAttention = isAttention  # biLSTM中是否用attention机制
        self.filter_sizes = filter_sizes  # 用于textCNN 卷积核大小
        self.num_filters = num_filters  # 用于textCNN 卷积核数量
        self.output_size = 128  # 用于RCNN模型中biLSTM 的输出维度大小
        self.num_blocks = num_blocks  # transformer encoder block的个数
        self.head = head  # multi-attention 的头数

        # 定义模型输入输出
        self.input_data = tf.placeholder(tf.int32, shape=[None, self.max_seq_len], name="input_data")
        self.y = tf.placeholder(tf.float32, shape=[None, self.num_classes], name="target")
        with tf.name_scope("embedding"):
            embedding = tf.Variable(tf.random_normal(shape=[self.vocab_size, self.embedding_size], dtype=tf.float32),
                                    name="word_embedding")
            if preEmbedding:
                embedding = embedding.assign(tf.cast(wordVector, tf.float32))
                tf.logging.info("已使用预训练的词向量")
        self.x = tf.nn.embedding_lookup(embedding, self.input_data)
        self.keep_prob = tf.placeholder(name="keep_prob", dtype=tf.float32)
        # print(self.x.shape)  # (?, 300, 50)
        with tf.name_scope("build_net"):
            self.build_net()

    def build_net(self):
        pass

    def get_metrics(self, y_pred, y_true):
        y_pred = np.argmax(y_pred, 1)
        y_true = np.argmax(y_true, 1)
        p = precision_score(y_true, y_pred)
        r = recall_score(y_true, y_pred)
        f = f1_score(y_true, y_pred)
        return p, r, f

class Dynamic_LSTM(TextClassifier):

    def build_net(self):
        # 定义softmax层参数
        weight = tf.get_variable("weight", shape=[self.n_hidden, self.num_classes], dtype=tf.float32,
                                 initializer=tf.random_normal_initializer)
        bias = tf.get_variable("bias", shape=[self.num_classes], dtype=tf.float32, initializer=tf.zeros_initializer)

        # 定义lstm层
        lstmcell = tf.nn.rnn_cell.BasicLSTMCell(self.n_hidden, forget_bias=1.0, state_is_tuple=True)
        lstmcell = tf.nn.rnn_cell.DropoutWrapper(lstmcell, output_keep_prob=self.keep_prob)
        lstmcell = tf.nn.rnn_cell.MultiRNNCell([lstmcell] * self.num_layers)
        outputs, states = tf.nn.dynamic_rnn(lstmcell, self.x, dtype=tf.float32)
        outputs = tf.transpose(outputs, [1, 0, 2])
        output = outputs[-1]
        output = tf.nn.dropout(output, keep_prob=self.keep_prob)
        self.pred = tf.matmul(output, weight) + bias
        self.loss = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=self.pred, labels=self.y))
        optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
        self.train_op = optimizer.minimize(self.loss)

class BiLSTMandAttention(TextClassifier):
    """双向LSTM+Attention"""

    def build_net(self):
        # self.x=tf.unstack(self.x,self.max_seq_len,1)
        tf.logging.info("the shape of input tensor:%s" % str(self.x.shape))
        # 定义softmax层参数

        with tf.name_scope("bilstm"):
            def cell_fw():
                cell = tf.nn.rnn_cell.BasicLSTMCell(self.n_hidden, forget_bias=1.0, state_is_tuple=True)
                return tf.nn.rnn_cell.DropoutWrapper(cell, output_keep_prob=self.keep_prob)

            def cell_bw():
                cell = tf.nn.rnn_cell.BasicLSTMCell(self.n_hidden, forget_bias=1.0, state_is_tuple=True)
                return tf.nn.rnn_cell.DropoutWrapper(cell, output_keep_prob=self.keep_prob)

            stack_fw_cell = tf.nn.rnn_cell.MultiRNNCell([cell_fw() for _ in range(self.num_layers)],
                                                        state_is_tuple=True)
            stack_bw_cell = tf.nn.rnn_cell.MultiRNNCell([cell_bw() for _ in range(self.num_layers)],
                                                        state_is_tuple=True)
            # outputs是正向和反向输出的元组,states是正向和反向最后时刻状态的元组
            outputs, states = tf.nn.bidirectional_dynamic_rnn(stack_fw_cell, stack_bw_cell, self.x, dtype=tf.float32)
        # 是否添加attention机制
        if self.isAttention:
            output = outputs[0] + outputs[1]
            with tf.name_scope("attention"):
                output = self.attention(output)  # 将各个时刻的状态值加权求和，结果做softmax分类
            weight = tf.get_variable("weight", shape=[self.n_hidden, self.num_classes], dtype=tf.float32,
                                     initializer=tf.random_normal_initializer)
            bias = tf.get_variable("bias", shape=[self.num_classes], dtype=tf.float32, initializer=tf.zeros_initializer)
        else:
            output = tf.concat([outputs[0], outputs[1]], axis=2)
            tf.logging.info("output:%s" % str(output.shape))
            output = tf.transpose(output, [1, 0, 2])
            output = output[-1]  # 将最后一个时间步做softmax分类
            tf.logging.info("output:%s" % str(output.shape))
            weight = tf.get_variable("weight", shape=[self.n_hidden * 2, self.num_classes], dtype=tf.float32,
                                     initializer=tf.random_normal_initializer)
            bias = tf.get_variable("bias", shape=[self.num_classes], dtype=tf.float32, initializer=tf.zeros_initializer)

        with tf.name_scope("outputs"):
            self.pred = tf.matmul(output, weight) + bias

        with tf.name_scope("loss_train_op"):
            self.loss = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=self.pred, labels=self.y))
            optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            self.train_op = optimizer.minimize(self.loss)

    def attention(self, H):
        # 初始化一个可训练权重向量，
        w = tf.Variable(tf.random_normal([self.n_hidden], stddev=0.1), name="attention_w")
        # 对输出做非线性转换
        M = tf.tanh(H)
        # 对w和M做矩阵运算，w=[batch_size,time_step,n_hidden],维度转换[batch_size*time_step,n_hidden]
        # 每一个时间步的输出向量转换为一个数字
        newM = tf.matmul(tf.reshape(M, [-1, self.n_hidden]), tf.reshape(w, [-1, 1]))
        # 对newM维度转换成[batch_size,time_step]
        restoreM = tf.reshape(newM, [-1, self.max_seq_len])
        # softmax归一化处理[batch_size,time_step]
        self.alpha = tf.nn.softmax(restoreM)

        # 利用求得的alpha的值对H进行加权求和，用矩阵运算直接操作
        H = tf.transpose(H, [0, 2, 1])
        r = tf.matmul(H, tf.reshape(self.alpha, [-1, self.max_seq_len, 1]))
        # 将三维压缩成二维sequeeze
        sequeezeR = tf.squeeze(r)
        output = tf.tanh(sequeezeR)
        # 做dropout处理
        output = tf.nn.dropout(output, keep_prob=self.keep_prob)
        return output

class TextCNN(TextClassifier):
    """TextCNN模型"""

    def build_net(self):
        # 扩展输入的维度
        self.x = tf.expand_dims(self.x, axis=-1, name="expand_dim")
        pooling_res = []
        for filter_size in self.filter_sizes:
            filter_shape = [filter_size, self.embedding_size, 1, self.num_filters]
            weight = tf.Variable(tf.random_normal(shape=filter_shape, stddev=1.0), dtype=tf.float32)
            bias = tf.Variable(tf.random_normal(shape=[self.num_filters], stddev=1.0), dtype=tf.float32)
            conv = tf.nn.conv2d(self.x, weight, [1, 1, 1, 1], "VALID", name="conv")
            relu = tf.nn.relu(tf.nn.bias_add(conv, bias), name="relu")
            max_pool = tf.nn.max_pool(relu, ksize=[1, self.max_seq_len - filter_size + 1, 1, 1], strides=[1, 1, 1, 1],
                                      padding="VALID", name="maxpool")
            pooling_res.append(max_pool)
        # 输入全连接层的神经元的个数
        num_fc = self.num_filters * len(self.filter_sizes)
        # print(pooling_res[0].shape)  # max_pool的shape (?, 1, 1, 256)
        fc_input = tf.concat(pooling_res, axis=3)  # 按照最后一个维度合并
        # print(fc_input.shape)  # 合并之后的shape (?, 1, 1, 1024)
        fc_input = tf.reshape(fc_input, shape=[-1, num_fc])
        # print(fc_input.shape)  # reshape之后的shape (?, 1024)
        fc_input = tf.nn.dropout(fc_input, keep_prob=self.keep_prob)
        with tf.name_scope("output"):
            w = tf.Variable(tf.truncated_normal(shape=[num_fc, self.num_classes], stddev=1.0), dtype=tf.float32,
                            name="w")
            b = tf.Variable(tf.zeros(shape=[self.num_classes]), name='b')
            self.pred = tf.matmul(fc_input, w) + b
        with tf.name_scope("loss_train"):
            self.loss = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=self.pred, labels=self.y))
            optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            self.train_op = optimizer.minimize(self.loss)

class RCNN(TextClassifier):
    """
    RNN+CNN 模型
    1.将双向LSTM的正向输出和反向输出和词向量拼成一个上下文向量，
    2.经tanh函数激活，送入最大池化层，
    3.然后将池化结果传入全连接做分类。
    """

    def build_net(self):
        with tf.name_scope("bilstm"):
            def cell_fw():
                cell = tf.nn.rnn_cell.BasicLSTMCell(self.n_hidden, forget_bias=1.0, state_is_tuple=True)
                return tf.nn.rnn_cell.DropoutWrapper(cell, output_keep_prob=self.keep_prob)

            def cell_bw():
                cell = tf.nn.rnn_cell.BasicLSTMCell(self.n_hidden, forget_bias=1.0, state_is_tuple=True)
                return tf.nn.rnn_cell.DropoutWrapper(cell, output_keep_prob=self.keep_prob)

            stack_fw_cell = tf.nn.rnn_cell.MultiRNNCell([cell_fw() for _ in range(self.num_layers)],
                                                        state_is_tuple=True)
            stack_bw_cell = tf.nn.rnn_cell.MultiRNNCell([cell_bw() for _ in range(self.num_layers)],
                                                        state_is_tuple=True)
            (output_fw, output_bw), states = tf.nn.bidirectional_dynamic_rnn(stack_fw_cell, stack_bw_cell, self.x,
                                                                             dtype=tf.float32)

        with tf.name_scope("context"):
            tf.logging.info("the shape of foreward vector:%s" % str(output_fw.shape))  # (?, 300, 128)
            shape = [tf.shape(output_fw)[0], 1, tf.shape(output_fw)[2]]
            # 获得上文信息 (?, 300, 128)[-1,max_seq_len,n_hidden]
            context_left = tf.concat([tf.zeros(shape), output_fw[:, :-1]], axis=1, name="contextleft")
            # 获得下文信息 (?, 300, 128)[-1,max_seq_len,n_hidden]
            context_right = tf.concat([output_bw[:, 1:], tf.zeros(shape)], axis=1, name="contextright")
            # 将前向/后向和当前词向量拼接在一起，作为最终的表征向量 (?, 300, 306)[-1,max_seq_len,2*n_hidden+embedding_size]
            context = tf.concat([context_left, self.x, context_right], axis=2, name="context")
            tf.logging.info("the shape of context_left:%s" % str(context_left.shape))
            tf.logging.info("the shape of context_right:%s" % str(context_right.shape))
            tf.logging.info("the shape of context%s" % str(context.shape))
        with tf.name_scope("rnn_output"):
            size = self.n_hidden * 2 + self.embedding_size
            weight = tf.Variable(tf.random_normal(shape=[size, self.output_size], stddev=1.0), dtype=tf.float32,
                                 name="rnn_output_weight")
            bias = tf.Variable(tf.random_normal(shape=[self.output_size]), dtype=tf.float32, name="rnn_output_bias")
            rnn_output = tf.tanh(tf.einsum("aij,jk->aik", context, weight) + bias)  # tf.einsum()方法用于指定维度消除
            tf.logging.info("the shape of rnn_output:%s" % str(rnn_output.shape))  # (?, 300, 128)
        with tf.name_scope("max_pool"):
            max_pool = tf.reduce_max(rnn_output, axis=1, name="max_pooling")
            tf.logging.info("the shape of max_pool output:%s" % str(max_pool.shape))  # (?, 128)将时间步的维度消除，取每个时间步最大的值
        with tf.name_scope("softmax"):
            w = tf.Variable(tf.truncated_normal([self.output_size, self.num_classes], stddev=1.0), dtype=tf.float32,
                            name="w")
            b = tf.Variable(tf.truncated_normal([self.num_classes]), dtype=tf.float32, name="b")

            pred = tf.matmul(max_pool, w) + b
            self.pred = tf.nn.softmax(pred)
        with tf.name_scope("loss_train"):
            l2_loss = tf.constant(0.0)
            l2_loss += tf.nn.l2_loss(w)
            l2_loss += tf.nn.l2_loss(b)
            # globalStep = tf.Variable(0, name="globalStep", trainable=False)
            self.loss = tf.reduce_mean(
                tf.nn.softmax_cross_entropy_with_logits(logits=self.pred, labels=self.y)) + 0.5 * l2_loss  # 增加L2正则
            optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            # gradAndtvar = optimizer.compute_gradients(self.loss)  # 获得梯度和变量
            # self.train_op = optimizer.apply_gradients(gradAndtvar, global_step=globalStep)
            self.train_op = optimizer.minimize(self.loss)

class CNN_GRU(TextClassifier):
    """CNN+GRU网络模型 将两者提取的信息concat，进行分类"""
    def build_net(self):
        with tf.name_scope("cnn"):
            x = tf.expand_dims(self.x, axis=-1, name="expand_dim")
            pooling_res = []
            for filter_size in self.filter_sizes:
                filter_shape = [filter_size, self.embedding_size, 1, self.num_filters]
                weight = tf.Variable(tf.truncated_normal(shape=filter_shape), dtype=tf.float32, name="weight")
                bias = tf.Variable(tf.zeros(shape=[self.num_filters]), dtype=tf.float32, name="bias")
                conv = tf.nn.conv2d(x, filter=weight, strides=[1, 1, 1, 1], padding="VALID", name="conv")
                relu = tf.nn.relu(tf.nn.bias_add(conv, bias))
                max_pool = tf.nn.max_pool(relu, ksize=[1, self.max_seq_len - filter_size + 1, 1, 1],
                                          strides=[1, 1, 1, 1], padding="VALID", name="pooling")
                pooling_res.append(max_pool)
            num_fc = self.num_filters * len(self.filter_sizes)
            fc_input = tf.concat(pooling_res, axis=3)
            cnn_output = tf.reshape(fc_input, shape=[-1, num_fc])
        with tf.name_scope("GRU"):
            cell = tf.nn.rnn_cell.GRUCell(self.n_hidden)
            cell = tf.nn.rnn_cell.DropoutWrapper(cell, output_keep_prob=self.keep_prob)
            cell = tf.nn.rnn_cell.MultiRNNCell([cell] * self.num_layers)
            outputs, _ = tf.nn.dynamic_rnn(cell, self.x, dtype=tf.float32)
            outputs = tf.transpose(outputs, [1, 0, 2])
            rnn_output = outputs[-1]
        with tf.name_scope("softmax"):
            fc = tf.concat([cnn_output, rnn_output], 1)
            fc = tf.nn.dropout(fc, keep_prob=self.keep_prob)
            w = tf.Variable(tf.truncated_normal(shape=[self.n_hidden + num_fc, self.num_classes]), dtype=tf.float32,
                            name="w")
            b = tf.Variable(tf.truncated_normal(shape=[self.num_classes]), dtype=tf.float32, name="b")
            self.pred = tf.matmul(fc, w) + b
        with tf.name_scope("loss_train"):
            l2_loss = tf.constant(0.0)
            l2_loss += tf.nn.l2_loss(w)
            l2_loss += tf.nn.l2_loss(b)
            self.loss = tf.reduce_mean(
                tf.nn.softmax_cross_entropy_with_logits(logits=self.pred, labels=self.y)) + 0.5 * l2_loss  # 增加L2正则
            optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            self.train_op = optimizer.minimize(self.loss)

class Transformer(TextClassifier):
    """
    transfromer encoder 用于文本分类
    1.position embedding，用one-hot，数据量小
    2.block层数设置为2
    3.sublayer加dropout正则化
    4.前馈层用卷积实现（全连接实现）待定
    """

    def _onehot_position_embedding(self, batch_size=64, max_len=300):
        """采取one-hot编码生成位置向量"""
        position_embedding = []
        for batch in range(batch_size):
            x = []
            for step in range(max_len):
                vector = np.zeros(max_len)
                vector[step] = 1
                x.append(vector)
            position_embedding.append(x)
        return np.array(position_embedding, dtype=np.float32)

    def _position_embedding(self):
        """采用论文中sin/cos生成位置向量"""
        batch_size = tf.shape(self.x)[0]
        # 生成位置索引，并扩张到所有的样本上
        positionIndex = tf.tile(tf.expand_dims(tf.range(self.max_seq_len), 0), [batch_size, 1])
        # 根据正弦和余弦函数获得每个位置上的第一部分
        positionEmbedding = np.array(
            [[pos / np.power(10000, (i - i % 2) / self.embedding_size) for i in range(self.embedding_size)] for pos in
             range(self.max_seq_len)])
        # 根据奇偶性分别用sin和cos函数包装
        positionEmbedding[:, 0::2] = np.sin(positionEmbedding[:, 0::2])
        positionEmbedding[:, 1::2] = np.cos(positionEmbedding[:, 1::2])
        # 转换成tensor格式
        positionEmbedding_ = tf.cast(positionEmbedding, dtype=tf.float32)
        # 获得三维的矩阵
        positionEmbedding = tf.nn.embedding_lookup(positionEmbedding_, positionIndex)
        return positionEmbedding

    def _multihead_attention(self,raw_inputs,queries,keys,size=None,causality=False):
        """
        参数
        :param raw_inputs:  原始输入，用于计算mask
        :param queries: 用于计算相关度，与key相同则是同一句子计算相关性
        :param keys: 用于计算相关度，word_embedding+position_embedding
        :param size: dense全连接层输出维度大小
        :param causality:
        :return:
        """
        if size is None:
            size=queries.shape[-1]
        #对数据分割映射(论文中是先分割在映射)
        Q=tf.layers.dense(queries,size,activation=tf.nn.relu)
        K=tf.layers.dense(keys,size,activation=tf.nn.relu)
        V=tf.layers.dense(keys,size,activation=tf.nn.relu)
        Q_=tf.concat(tf.split(Q,self.head,axis=-1),axis=0)
        K_=tf.concat(tf.split(K,self.head,axis=-1),axis=0)
        V_=tf.concat(tf.split(V,self.head,axis=-1),axis=0)
        #输出维度
        tf.logging.info("the shape of Q:%s"%str(Q.shape)) #(64, 300, 350)
        tf.logging.info("the shape of Q_:%s"%str(Q_.shape)) #(320, 300, 70)

        #计算k和q的点积，维度[batch_size*numheads,q_len,k_len]
        similary=tf.matmul(Q_,tf.transpose(K_,[0,2,1]))
        #进行缩放处理
        scaled_similary=similary/(K_.get_shape().as_list()[-1]**0.5)
        #对填充词的处理，当padding为0是时，计算的权重应该是0，需要将其mask为0，则需要q或者k一方为0
        #将每一个时序上的向量中的值相加取平均值
        key_mask=tf.sign(tf.abs(tf.reduce_sum(raw_inputs,axis=-1))) #[batch_size,time_step]
        #利用tf.tile进行张量扩张，维度[batch_size*numhead,key_len]
        key_mask=tf.tile(key_mask,[self.head,1])
        #增加一个维度，并进行扩张，维度[batch_size*numhead,q_len,k_len]
        key_mask=tf.tile(tf.expand_dims(key_mask,1),[1,queries.shape[1],1])
        tf.logging.info("the shape of key_mask:%s"%str(key_mask.shape)) #(?, 300, 300)

        #生成与scale_similary相同维度的tensor，然后得到负无穷大的值
        paddings=tf.ones_like(scaled_similary)*(-2**(32+1))
        #如果keymask中的值为0则用padding替换，
        masked_similary=tf.where(tf.equal(key_mask,0),paddings,scaled_similary)

        #计算权重系数
        weight=tf.nn.softmax(masked_similary)
        #加权求和，维度[batch_size*numheads,time_step,embedding_size/numheads]
        outputs=tf.matmul(weight,V_)
        #将outputs重组
        outputs=tf.concat(tf.split(outputs,self.head,0),axis=2)
        outputs=tf.nn.dropout(outputs,keep_prob=self.keep_prob)

        #对每一个sublayer建立残差连接
        outputs+=queries
        #标准层
        outputs = self._layer_normalization(outputs)
        return outputs

    def _layer_normalization(self,inputs):
        """Batch Normalization层"""
        inputshape=inputs.shape
        paramshape=inputshape[-1]
        #LayerNorm是在最后的维度上计算输入的数据的均值和方差，维度[batch_size,time_step,1]
        mean,var=tf.nn.moments(inputs,[-1],keep_dims=True)

        beta=tf.Variable(tf.zeros(paramshape))
        gamma=tf.Variable(tf.ones(paramshape))
        normalized=(inputs-mean)/((var+1.0)**0.5)
        outputs=gamma*normalized+beta   #重构
        return outputs

    def _feed_forward(self,inputs):
        filters=[128,self.embedding_size+self.max_seq_len] #[内层，外层]
        #前馈层采用卷积神经网络 1维卷积 实际的维度还是二维的[batch_size,time_step,embedding_size]
        param={"inputs":inputs,"filters":filters[0],"kernel_size":1,"activation":tf.nn.relu,"use_bias":True}
        outputs=tf.layers.conv1d(**param)
        param={"inputs":outputs,"filters": filters[1], "kernel_size": 1,"activation": None, "use_bias": True}
        outputs=tf.layers.conv1d(**param)

        #残差连接和归一化处理
        outputs+=inputs
        outputs=self._layer_normalization(outputs)
        return outputs

    def build_net(self):
        with tf.name_scope("input"):
            position_embedding = self._onehot_position_embedding(batch_size=64, max_len=self.max_seq_len)
            embedded_words = tf.concat([self.x, position_embedding], axis=-1)
            tf.logging.info("the shape of position embedding:%s" % str(position_embedding.shape))
            tf.logging.info("the shape of embedded words:%s" % str(embedded_words.shape))

        with tf.name_scope("encoder"):
            for i in range(self.num_blocks):
                with tf.name_scope("encoder-%d"%(i+1)):
                    multihead_attention=self._multihead_attention(raw_inputs=self.x,queries=embedded_words,keys=embedded_words)
                    embedded_words=self._feed_forward(inputs=multihead_attention)
            tf.logging.info("the last encoder block's shape:%s"%str(embedded_words.shape)) #(64, 300, 350)

            outputs=tf.reshape(embedded_words,[-1,self.max_seq_len*(self.embedding_size+self.max_seq_len)])
            output_size=outputs.shape[-1]

        with tf.name_scope("dropout"):
            outputs=tf.nn.dropout(outputs,keep_prob=self.keep_prob)

        with tf.name_scope("output"):
            w=tf.get_variable("output_w",shape=[output_size,self.num_classes],initializer=tf.random_normal_initializer)
            b=tf.get_variable("output_b",shape=[self.num_classes],initializer=tf.zeros_initializer)
            self.pred=tf.matmul(outputs,w)+b

        with tf.name_scope("loss_train"):
            self.loss=tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=self.pred,labels=self.y))
            optimizer=tf.train.AdamOptimizer(learning_rate=self.learning_rate)
            self.train_op=optimizer.minimize(self.loss)

class Capsule(object):
    """
    Capsule应用于文本分类，
    论文《Investigating Capsule Networks with Dynamic Routing for Text Classification》
    第一层：卷积层  卷积核K1=3，卷积核个数B=32
    第二层：基础capsule  卷积核个数C=32，输出g(W*Mi+b)
    """
    def __init__(self,
                 max_seq_len=300,
                 embedding_size=50,
                 num_classes=2,
                 vocab_size=None,
                 wordVector=None,
                 preEmbedding=True,
                 filter_size=3,
                 num_filters=32,
                 learning_rate=0.0005,
                 epoch=20):
        self.max_seq_len=max_seq_len  #最大序列长度
        self.embedding_size=embedding_size  #词向量维度大小
        self.num_classes=num_classes  #类别标签个数
        self.vocab_size=vocab_size  #embedding中词的个数
        self.filter_size=filter_size  #卷积核大小，论文中大小为3
        self.num_filters=num_filters  #卷积核数量
        self.lr=learning_rate
        self.epoch=epoch

        self.keep_prob=tf.placeholder(tf.float32,name="keep_prob")
        self.input_data = tf.placeholder(tf.int32, shape=[None, self.max_seq_len], name="input_data")
        self.y = tf.placeholder(tf.float32, shape=[None, self.num_classes], name="target")
        with tf.name_scope("embedding"):
            embedding = tf.Variable(tf.random_normal(shape=[self.vocab_size, self.embedding_size], dtype=tf.float32),
                                    name="word_embedding")
            if preEmbedding:
                embedding = embedding.assign(tf.cast(wordVector, tf.float32))
                tf.logging.info("已使用预训练的词向量")
        x = tf.nn.embedding_lookup(embedding, self.input_data)
        self.x=tf.expand_dims(x,-1)

        with tf.name_scope("build_net"):
            self.build_net()

    def build_net(self):
        with tf.name_scope("conv_layer"):
            filter_shape=[self.filter_size,self.embedding_size,1,self.num_filters]  #[3,50,1,32]
            w_conv=tf.Variable(tf.truncated_normal(shape=filter_shape),dtype=tf.float32,name="w_conv")
            b_conv=tf.Variable(tf.truncated_normal(shape=[self.num_filters]),dtype=tf.float32,name="b_conv")
            conv=tf.nn.conv2d(self.x,filter=w_conv,strides=[1,2,1,1],padding="VALID",name="conv")
            conv_output=tf.nn.relu(tf.nn.bias_add(conv,b_conv))
            #返回shape[batch_size,max_len-filter_size+1,1,num_filters]
            tf.logging.info("the shape of conv_output {}".format(conv_output.get_shape()))

        with tf.name_scope("primary_capsule_layer"):
            w_pc=tf.Variable(tf.truncated_normal(shape=[1,1,32,16]),dtype=tf.float32,name="w_pc")
            b_pc=tf.Variable(tf.truncated_normal(shape=[16]),dtype=tf.float32,name="b_pc")
            conv=tf.nn.conv2d(conv_output,filter=w_pc,strides=[1,1,1,1],padding="VALID",name="pc_conv")
            pose=tf.nn.relu(tf.nn.bias_add(conv,b_pc),name="pc_relu")
            pose_shape=pose.get_shape().as_list()
            pose=tf.reshape(pose,[-1, pose_shape[1], pose_shape[2], 16, 16])
            tf.logging.info("the shape of primary_pose {}".format(pose.get_shape())) #(?, 149, 1, 16, 16)
            pose=self.squash(pose)
            bias_2=tf.Variable(tf.truncated_normal(shape=[1,16]),dtype=tf.float32,name="bias_2")
            activation=tf.sqrt(np.sum(np.square(pose),axis=-1))+bias_2

        with tf.name_scope("capsule_conv_layer"):
            inputs_pose_shape=pose.get_shape().as_list()
            shape=[3,1,16,16]
            hk_offsets=[
                [int(h_offset+k_offset) for k_offset in range(0,shape[0])] for h_offset in range(0,inputs_pose_shape[1]+1-shape[0],1)
            ]
            wk_offsets = [
                [int(w_offset + k_offset) for k_offset in range(0, shape[1])] for w_offset in range(0, inputs_pose_shape[2] + 1 - shape[1], 1)
            ]
            print(hk_offsets)
            print(wk_offsets)
            inputs_poses_patches = tf.transpose(
                tf.gather(
                    tf.gather(
                        pose, hk_offsets, axis=1, name='gather_poses_height_kernel'
                    ), wk_offsets, axis=3, name='gather_poses_width_kernel'
                ), perm=[0, 1, 3, 2, 4, 5, 6], name='inputs_poses_patches'
            )
            tf.logging.info('i_poses_patches shape: {}'.format(inputs_poses_patches.get_shape()))

            inputs_poses_shape = inputs_poses_patches.get_shape().as_list()
            inputs_poses_patches = tf.reshape(inputs_poses_patches, [
                -1, shape[0] * shape[1] * shape[2], inputs_poses_shape[-1]
            ])
            i_activations_patches = tf.transpose(
                tf.gather(
                    tf.gather(
                        activation, hk_offsets, axis=1, name='gather_activations_height_kernel'
                    ), wk_offsets, axis=3, name='gather_activations_width_kernel'
                ), perm=[0, 1, 3, 2, 4, 5, 6], name='inputs_activations_patches'
            )
            tf.logging.info('i_activations_patches shape: {}'.format(i_activations_patches.get_shape()))
            i_activations_patches = tf.reshape(i_activations_patches, [-1, shape[0] * shape[1] * shape[2]])
            u_hat_vecs = self.vec_transformationByConv(inputs_poses_patches,
                                                       inputs_poses_shape[-1],
                                                       shape[0] * shape[1] * shape[2],
                                                       inputs_poses_shape[-1],
                                                       shape[3])
            tf.logging.info('capsule conv votes shape: {}'.format(u_hat_vecs.get_shape()))
            beta_a=tf.Variable(tf.random_normal(shape=[1,shape[3]]),dtype=tf.float32,name="beta_a")
            poses, activations = self.routing(u_hat_vecs, beta_a, 3, shape[3], i_activations_patches)
            poses = tf.reshape(poses, [
                -1, inputs_poses_shape[1],
                inputs_poses_shape[2], shape[3],
                inputs_poses_shape[-1]]
                               )
            activations = tf.reshape(activations, [
                -1, inputs_poses_shape[1],
                inputs_poses_shape[2], shape[3]])
            tf.logging.info("capsule conv poses dimension:{}".format(poses.get_shape()))
            tf.logging.info("capsule conv activations dimension:{}".format(activations.get_shape()))

        with tf.name_scope("capsule_flatten"):
            input_pose_shape=poses.get_shape().as_list()
            poses = tf.reshape(poses, [
                -1, input_pose_shape[1] * input_pose_shape[2] * input_pose_shape[3], input_pose_shape[-1]])
            activations = tf.reshape(activations, [
                -1, input_pose_shape[1] * input_pose_shape[2] * input_pose_shape[3]])
            tf.logging.info("flatten poses dimension:{}".format(poses.get_shape()))
            tf.logging.info("flatten activations dimension:{}".format(activations.get_shape()))

        with tf.name_scope("capsule_fc_layer"):
            input_pose_shape=poses.get_shape().as_list()
            u_hat_vecs = self.vec_transformationByConv(
                poses,
                input_pose_shape[-1], input_pose_shape[1],
                input_pose_shape[-1], self.num_classes,
            )
            tf.logging.info('votes shape: {}'.format(u_hat_vecs.get_shape()))
            beta_b=tf.Variable(tf.random_normal(shape=[1,self.num_classes]),tf.float32,name="beta_b")
            poses, activations = self.routing(u_hat_vecs, beta_b, 3, self.num_classes,i_activations=None)
            tf.logging.info('capsule fc shape: {}'.format(poses.get_shape()))
            # print(activations)
            # exit()
        with tf.name_scope("loss_train"):
            self.pred=activations
            self.loss=tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=self.pred,labels=self.y))
            optimizer=tf.train.AdamOptimizer(learning_rate=self.lr)
            gradAndVar=optimizer.compute_gradients(self.loss)
            self.train_op=optimizer.apply_gradients(grads_and_vars=gradAndVar)

    def squash(self,x):
        """向量归一化，使长度小于1"""
        s_squared_norm = np.sum(np.square(x), axis=-1, keepdims=True)
        scale =tf.sqrt(s_squared_norm) / (0.5 + s_squared_norm)
        return scale * x

    def vec_transformationByConv(self,poses, input_capsule_dim, input_capsule_num, output_capsule_dim, output_capsule_num):
        kernel=tf.Variable(tf.random_normal(shape=[1, input_capsule_dim, output_capsule_dim * output_capsule_num]),dtype=tf.float32,name="weight")
        tf.logging.info('poses: {}'.format(poses.get_shape()))
        tf.logging.info('kernel: {}'.format(kernel.get_shape()))
        u_hat_vecs = tf.nn.conv1d(poses, kernel,stride=1,padding="VALID")
        u_hat_vecs = tf.reshape(u_hat_vecs, (-1, input_capsule_num, output_capsule_num, output_capsule_dim))
        u_hat_vecs = tf.transpose(u_hat_vecs, (0, 2, 1, 3))
        return u_hat_vecs

    def routing(self, u_hat_vecs, beta_a, iterations, output_capsule_num, i_activations):
        from keras import backend as K
        b = tf.zeros_like(u_hat_vecs[:, :, :, 0])
        if i_activations is not None:
            i_activations = i_activations[..., tf.newaxis]
        for i in range(iterations):
            if False:
                leak = tf.zeros_like(b, optimize=True)
                leak = tf.reduce_sum(leak, axis=1, keep_dims=True)
                leaky_logits = tf.concat([leak, b], axis=1)
                leaky_routing = tf.nn.softmax(leaky_logits, dim=1)
                c = tf.split(leaky_routing, [1, output_capsule_num], axis=1)[1]
            else:
                c = tf.nn.softmax(b, 1)
                #        if i_activations is not None:
            #            tf.transpose(tf.transpose(c, perm=[0,2,1]) * i_activations, perm=[0,2,1])
            outputs = self.squash(K.batch_dot(c, u_hat_vecs, [2, 2]))
            if i < iterations - 1:
                b = b + K.batch_dot(outputs, u_hat_vecs, [2, 3])
        poses = outputs
        activations = K.sqrt(K.sum(K.square(poses), 2))
        return poses, activations

max_stars_repo_name: longshuicui/IMDB-movie-reviews-analysis

max_issues_repo_name: longshuicui/IMDB-movie-reviews-analysis

max_forks_repo_name: longshuicui/IMDB-movie-reviews-analysis


text:
import numpy as np
from scipy.fftpack import ifftshift, ifft2, fft2
from scipy.interpolate import RegularGridInterpolator

from aspyre.utils import ensure
from aspyre.utils.math import grid_2d
from aspyre.utils.fft import centered_fft2, centered_ifft2
from aspyre.utils.matrix import roll_dim, unroll_dim

# TODO: The following functions all need to move inside the Image class
def im_translate(im, shifts):
    """
    Translate image by shifts
    :param im: An array of size L-by-L-by-n containing images to be translated.
    :param shifts: An array of size 2-by-n specifying the shifts in pixels.
        Alternatively, it can be a column vector of length 2, in which case the same shifts is applied to each image.
    :return: The images translated by the shifts, with periodic boundaries.
    """

    n_im = im.shape[-1]
    n_shifts = shifts.shape[-1]

    ensure(shifts.shape[0] == 2, "shifts must be 2xn")
    ensure(n_shifts == 1 or n_shifts == n_im, "no. of shifts must be 1 or match the no. of images")
    ensure(im.shape[0] == im.shape[1], "images must be square")

    L = im.shape[0]
    im_f = fft2(im, axes=(0, 1))
    grid_1d = ifftshift(np.ceil(np.arange(-L/2, L/2))) * 2 * np.pi / L
    om_x, om_y = np.meshgrid(grid_1d, grid_1d, indexing='ij')

    phase_shifts_x = np.broadcast_to(-shifts[0, :], (L, L, n_shifts))
    phase_shifts_y = np.broadcast_to(-shifts[1, :], (L, L, n_shifts))
    phase_shifts = (om_x[:, :, np.newaxis] * phase_shifts_x) + (om_y[:, :, np.newaxis] * phase_shifts_y)

    mult_f = np.exp(-1j * phase_shifts)
    im_translated_f = im_f * mult_f
    im_translated = ifft2(im_translated_f, axes=(0, 1))
    im_translated = np.real(im_translated)

    return im_translated

def im_downsample(im, L_ds):
    """
    Blur and downsample image
    :param im: Set of images to be downsampled in the form of an array L-by-L-by-K, where K is the number of images.
    :param L_ds: The desired resolution of the downsampled images. Must be smaller than L.
    :return: An array of the form L_ds-by-L_ds-by-K consisting of the blurred and downsampled images.
    """
    N = im.shape[0]
    grid = grid_2d(N)
    grid_ds = grid_2d(L_ds)

    im_ds = np.zeros((L_ds, L_ds, im.shape[2])).astype(im.dtype)

    # x, y values corresponding to 'grid'. This is what scipy interpolator needs to function.
    x = y = np.ceil(np.arange(-N/2, N/2)) / (N/2)

    mask = (np.abs(grid['x']) < L_ds/N) & (np.abs(grid['y']) < L_ds/N)
    im = np.real(centered_ifft2(centered_fft2(im) * np.expand_dims(mask, 2)))

    for s in range(im_ds.shape[-1]):
        interpolator = RegularGridInterpolator(
            (x, y),
            im[:, :, s],
            bounds_error=False,
            fill_value=0
        )
        im_ds[:, :, s] = interpolator(np.dstack([grid_ds['x'], grid_ds['y']]))

    return im_ds

def im_filter(im, filt, *args, **kwargs):
    # TODO: Move inside appropriate object
    L = im.shape[0]
    im, sz_roll = unroll_dim(im, 3)
    filter_vals = filt.evaluate_grid(L, *args, **kwargs)
    im_f = centered_fft2(im)
    if im_f.ndim > filter_vals.ndim:
        im_f = np.expand_dims(filter_vals, 2) * im_f
    else:
        im_f = filter_vals * im_f
    im = centered_ifft2(im_f)
    im = np.real(im)
    im = roll_dim(im, sz_roll)

    return im

ext: py

max_stars_repo_name: ComputationalCryoEM/ASPIRE

max_issues_repo_name: ComputationalCryoEM/ASPIRE

max_forks_repo_name: computationalcryoem/aspyre

reason: import numpy,from scipy


text:
from collections import defaultdict
import math
import numpy as np

from krwordrank.graph import hits

def summarize_with_keywords(texts, num_keywords=100, stopwords=None, min_count=5,
    max_length=10, beta=0.85, max_iter=10, num_rset=-1, verbose=False):
    """
    It train KR-WordRank to extract keywords from texts.

        >>> from krwordrank.word import summarize_with_keywords

        >>> texts = [] # list of str
        >>> keywords = summarize_with_keywords(texts, num_keywords=100, min_count=5)

    Arguments
    ---------
    texts : list of str
        Each str is a sentence.
    num_keywords : int
        Number of keywords extracted from KR-WordRank
        Default is 100.
    stopwords : None or set of str
        Stopwords list for keyword and key-sentence extraction
    min_count : int
        Minimum frequency of subwords used to construct subword graph
        Default is 5
    max_length : int
        Maximum length of subwords used to construct subword graph
        Default is 10
    beta : float
        PageRank damping factor. 0 < beta < 1
        Default is 0.85
    max_iter : int
        Maximum number of iterations of HITS algorithm.
        Default is 10
    num_rset : int
        Number of R set words sorted by rank. It will be used to L-part word filtering.
        Default is -1.
    verbose : Boolean
        If True, it shows training status
        Default is False

    Returns
    -------
    keywords : dict
        Word : rank dictionary. keywords[str] = float

    Usage
    -----
        >>> from krwordrank.word import summarize_with_keywords

        >>> texts = [] # list of str
        >>> keywords = summarize_with_keywords(texts, num_keywords=100, min_count=5)
    """
    # train KR-WordRank
    wordrank_extractor = KRWordRank(
        min_count = min_count,
        max_length = max_length,
        verbose = verbose
        )

    keywords, rank, graph = wordrank_extractor.extract(texts,
        beta, max_iter, num_rset=num_rset)

    # stopword filtering
    if stopwords is None:
        stopwords = {}
    keywords = {word:r for word, r in keywords.items() if not (word in stopwords)}

    # top rank filtering
    if num_keywords > 0:
        keywords = {word:r for word, r in sorted(keywords.items(), key=lambda x:-x[1])[:num_keywords]}

    return keywords

class KRWordRank:
    """Unsupervised Korean Keyword Extractor

    Implementation of Kim, H. J., Cho, S., & Kang, P. (2014). KR-WordRank:
    An Unsupervised Korean Word Extraction Method Based on WordRank.
    Journal of Korean Institute of Industrial Engineers, 40(1), 18-33.

    Arguments
    ---------
    min_count : int
        Minimum frequency of subwords used to construct subword graph
        Default is 5
    max_length : int
        Maximum length of subwords used to construct subword graph
        Default is 10
    verbose : Boolean
        If True, it shows training status
        Default is False

    Usage
    -----
        >>> from krwordrank.word import KRWordRank

        >>> texts = ['예시 문장 입니다', '여러 문장의 list of str 입니다', ... ]
        >>> wordrank_extractor = KRWordRank()
        >>> keywords, rank, graph = wordrank_extractor.extract(texts, beta, max_iter, verbose)
    """
    def __init__(self, min_count=5, max_length=10, verbose=False):
        self.min_count = min_count
        self.max_length = max_length
        self.verbose = verbose
        self.sum_weight = 1
        self.vocabulary = {}
        self.index2vocab = []

    def scan_vocabs(self, docs):
        """
        It scans subwords positioned of left-side (L) and right-side (R) of words.
        After scanning was done, KR-WordRank has index2vocab as class attribute.

        Arguments
        ---------
        docs : list of str
            Sentence list

        Returns
        -------
        counter : dict
            {(subword, 'L')] : frequency}
        """
        self.vocabulary = {}
        if self.verbose:
            print('scan vocabs ... ')

        counter = {}
        for doc in docs:

            for token in doc.split():
                len_token = len(token)
                counter[(token, 'L')] = counter.get((token, 'L'), 0) + 1

                for e in range(1, min(len(token), self.max_length)):
                    if (len_token - e) > self.max_length:
                        continue

                    l_sub = (token[:e], 'L')
                    r_sub = (token[e:], 'R')
                    counter[l_sub] = counter.get(l_sub, 0) + 1
                    counter[r_sub] = counter.get(r_sub, 0) + 1

        counter = {token:freq for token, freq in counter.items() if freq >= self.min_count}
        for token, _ in sorted(counter.items(), key=lambda x:x[1], reverse=True):
            self.vocabulary[token] = len(self.vocabulary)

        self._build_index2vocab()

        if self.verbose:
            print('num vocabs = %d' % len(counter))
        return counter

    def _build_index2vocab(self):
        self.index2vocab = [vocab for vocab, index in sorted(self.vocabulary.items(), key=lambda x:x[1])]
        self.sum_weight = len(self.index2vocab)

    def extract(self, docs, beta=0.85, max_iter=10, num_keywords=-1,
        num_rset=-1, vocabulary=None, bias=None, rset=None):
        """
        It constructs word graph and trains ranks of each node using HITS algorithm.
        After training it selects suitable subwords as words.

        Arguments
        ---------
        docs : list of str
            Sentence list.
        beta : float
            PageRank damping factor. 0 < beta < 1
            Default is 0.85
        max_iter : int
            Maximum number of iterations of HITS algorithm.
            Default is 10
        num_keywords : int
            Number of keywords sorted by rank.
            Default is -1. If the vaule is negative, it returns all extracted words.
        num_rset : int
            Number of R set words sorted by rank. It will be used to L-part word filtering.
            Default is -1.
        vocabulary : None or dict
            User specified vocabulary to index mapper
        bias : None or dict
            User specified HITS bias term
        rset : None or dict
            User specfied R set

        Returns
        -------
        keywords : dict
            word : rank dictionary. {str:float}
        rank : dict
            subword : rank dictionary. {int:float}
        graph : dict of dict
            Adjacent subword graph. {int:{int:float}}

        Usage
        -----
            >>> from krwordrank.word import KRWordRank

            >>> texts = ['예시 문장 입니다', '여러 문장의 list of str 입니다', ... ]
            >>> wordrank_extractor = KRWordRank()
            >>> keywords, rank, graph = wordrank_extractor.extract(texts, beta, max_iter, verbose)
        """

        rank, graph = self.train(docs, beta, max_iter, vocabulary, bias)

        lset = {self.int2token(idx)[0]:r for idx, r in rank.items() if self.int2token(idx)[1] == 'L'}
        if not rset:
            rset = {self.int2token(idx)[0]:r for idx, r in rank.items() if self.int2token(idx)[1] == 'R'}

        if num_rset > 0:
            rset = {token:r for token, r in sorted(rset.items(), key=lambda x:-x[1])[:num_rset]}

        keywords = self._select_keywords(lset, rset)
        keywords = self._filter_compounds(keywords)
        keywords = self._filter_subtokens(keywords)

        if num_keywords > 0:
            keywords = {token:r for token, r in sorted(keywords.items(), key=lambda x:-x[1])[:num_keywords]}

        return keywords, rank, graph

    def _select_keywords(self, lset, rset):
        keywords = {}
        for word, r in sorted(lset.items(), key=lambda x:x[1], reverse=True):
            len_word = len(word)
            if len_word == 1:
                continue

            is_compound = False
            for e in range(2, len_word):
                if (word[:e] in keywords) and (word[:e] in rset):
                    is_compound = True
                    break

            if not is_compound:
                keywords[word] = r

        return keywords

    def _filter_compounds(self, keywords):
        keywords_= {}
        for word, r in sorted(keywords.items(), key=lambda x:x[1], reverse=True):
            len_word = len(word)

            if len_word <= 2:
                keywords_[word] = r
                continue

            if len_word == 3:
                if word[:2] in keywords_:
                    continue

            is_compound = False
            for e in range(2, len_word - 1):
                # fixed. comment from Y. cho
                if (word[:e] in keywords) and (word[e:] in keywords):
                    is_compound = True
                    break

            if not is_compound:
                keywords_[word] = r

        return keywords_

    def _filter_subtokens(self, keywords):
        subtokens = set()
        keywords_ = {}

        for word, r in sorted(keywords.items(), key=lambda x:x[1], reverse=True):
            subs = {word[:e] for e in range(2, len(word)+1)}

            is_subtoken = False
            for sub in subs:
                if sub in subtokens:
                    is_subtoken = True
                    break

            if not is_subtoken:
                keywords_[word] = r
                subtokens.update(subs)

        return keywords_

    def train(self, docs, beta=0.85, max_iter=10, vocabulary=None, bias=None):
        """
        It constructs word graph and trains ranks of each node using HITS algorithm.
        Use this function only when you want to train rank of subwords

        Arguments
        ---------
        docs : list of str
            Sentence list.
        beta : float
            PageRank damping factor. 0 < beta < 1
            Default is 0.85
        max_iter : int
            Maximum number of iterations of HITS algorithm.
            Default is 10
        vocabulary : None or dict
            User specified vocabulary to index mapper
        bias : None or dict
            User specified HITS bias term

        Returns
        -------
        rank : dict
            subword : rank dictionary. {int:float}
        graph : dict of dict
            Adjacent subword graph. {int:{int:float}}
        """
        if (not vocabulary) and (not self.vocabulary):
            self.scan_vocabs(docs)
        elif (not vocabulary):
            self.vocabulary = vocabulary
            self._build_index2vocab()

        graph = self._construct_word_graph(docs)

        rank = hits(graph, beta, max_iter, bias,
                    sum_weight=self.sum_weight,
                    number_of_nodes=len(self.vocabulary),
                    verbose=self.verbose
                    )

        return rank, graph

    def token2int(self, token):
        """
        Arguments
        ---------
        token : tuple
            (subword, 'L') or (subword, 'R')
            For example, ('이것', 'L') or ('은', 'R')

        Returns
        -------
        index : int
            Corresponding index
            If it is unknown, it returns -1
        """
        return self.vocabulary.get(token, -1)

    def int2token(self, index):
        """
        Arguments
        ---------
        index : int
            Token index

        Returns
        -------
        token : tuple
            Corresponding index formed such as (subword, 'L') or (subword, 'R')
            For example, ('이것', 'L') or ('은', 'R').
            If it is unknown, it returns None
        """
        return self.index2vocab[index] if (0 <= index < len(self.index2vocab)) else None

    def _construct_word_graph(self, docs):
        def normalize(graph):
            graph_ = defaultdict(lambda: defaultdict(lambda: 0))
            for from_, to_dict in graph.items():
                sum_ = sum(to_dict.values())
                for to_, w in to_dict.items():
                    graph_[to_][from_] = w / sum_
            graph_ = {t:dict(fd) for t, fd in graph_.items()}
            return graph_

        graph = defaultdict(lambda: defaultdict(lambda: 0))
        for doc in docs:

            tokens = doc.split()

            if not tokens:
                continue

            links = []
            for token in tokens:
                links += self._intra_link(token)

            if len(tokens) > 1:
                tokens = [tokens[-1]] + tokens + [tokens[0]]
                links += self._inter_link(tokens)

            links = self._check_token(links)
            if not links:
                continue

            links = self._encode_token(links)
            for l_node, r_node in links:
                graph[l_node][r_node] += 1
                graph[r_node][l_node] += 1

        # reverse for inbound graph. but it normalized with sum of outbound weight
        graph = normalize(graph)
        return graph

    def _intra_link(self, token):
        links = []
        len_token = len(token)
        for e in range(1, min(len_token, 10)):
            if (len_token - e) > self.max_length:
                continue
            links.append( ((token[:e], 'L'), (token[e:], 'R')) )
        return links

    def _inter_link(self, tokens):
        def rsub_to_token(t_left, t_curr):
            return [((t_left[-b:], 'R'), (t_curr, 'L')) for b in range(1, min(10, len(t_left)))]
        def token_to_lsub(t_curr, t_rigt):
            return [((t_curr, 'L'), (t_rigt[:e], 'L')) for e in range(1, min(10, len(t_rigt)))]

        links = []
        for i in range(1, len(tokens)-1):
            links += rsub_to_token(tokens[i-1], tokens[i])
            links += token_to_lsub(tokens[i], tokens[i+1])
        return links

    def _check_token(self, token_list):
        return [(token[0], token[1]) for token in token_list if (token[0] in self.vocabulary and token[1] in self.vocabulary)]

    def _encode_token(self, token_list):
        return [(self.vocabulary[token[0]],self.vocabulary[token[1]]) for token in token_list]

reason: import numpy


text:
# Copyright (c) 2019, NVIDIA CORPORATION. All rights reserved.
#
# Redistribution and use in source and binary forms, with or without
# modification, are permitted provided that the following conditions
# are met:
#  * Redistributions of source code must retain the above copyright
#    notice, this list of conditions and the following disclaimer.
#  * Redistributions in binary form must reproduce the above copyright
#    notice, this list of conditions and the following disclaimer in the
#    documentation and/or other materials provided with the distribution.
#  * Neither the name of NVIDIA CORPORATION nor the names of its
#    contributors may be used to endorse or promote products derived
#    from this software without specific prior written permission.
#
# THIS SOFTWARE IS PROVIDED BY THE COPYRIGHT HOLDERS ``AS IS'' AND ANY
# EXPRESS OR IMPLIED WARRANTIES, INCLUDING, BUT NOT LIMITED TO, THE
# IMPLIED WARRANTIES OF MERCHANTABILITY AND FITNESS FOR A PARTICULAR
# PURPOSE ARE DISCLAIMED.  IN NO EVENT SHALL THE COPYRIGHT OWNER OR
# CONTRIBUTORS BE LIABLE FOR ANY DIRECT, INDIRECT, INCIDENTAL, SPECIAL,
# EXEMPLARY, OR CONSEQUENTIAL DAMAGES (INCLUDING, BUT NOT LIMITED TO,
# PROCUREMENT OF SUBSTITUTE GOODS OR SERVICES; LOSS OF USE, DATA, OR
# PROFITS; OR BUSINESS INTERRUPTION) HOWEVER CAUSED AND ON ANY THEORY
# OF LIABILITY, WHETHER IN CONTRACT, STRICT LIABILITY, OR TORT
# (INCLUDING NEGLIGENCE OR OTHERWISE) ARISING IN ANY WAY OUT OF THE USE
# OF THIS SOFTWARE, EVEN IF ADVISED OF THE POSSIBILITY OF SUCH DAMAGE.

import sys
sys.path.append("../common")

import argparse
from builtins import range
from builtins import str
from future.utils import iteritems
import os
import time
import threading
import traceback
import numpy as np
import test_util as tu
from tensorrtserver.api import *
import tensorrtserver.api.server_status_pb2 as server_status

FLAGS = None
CORRELATION_ID_BLOCK_SIZE = 100
DEFAULT_TIMEOUT_MS = 5000
SEQUENCE_LENGTH_MEAN = 16
SEQUENCE_LENGTH_STDEV = 8

_thread_exceptions = []
_thread_exceptions_mutex = threading.Lock()

class TimeoutException(Exception):
    pass

def check_sequence_async(ctx, trial, model_name, input_dtype, steps,
                         timeout_ms=DEFAULT_TIMEOUT_MS, batch_size=1, sequence_name="<unknown>"):
    """Perform sequence of inferences using async run. The 'steps' holds
    a list of tuples, one for each inference with format:

    (flag_str, value, expected_result, delay_ms)

    """
    if (("savedmodel" in trial) or ("graphdef" in trial) or
        ("netdef" in trial) or ("custom" in trial) or
        ("plan" in trial)):
        tensor_shape = (1,)
    else:
        assert False, "unknown trial type: " + trial

    # Execute the sequence of inference...
    seq_start_ms = int(round(time.time() * 1000))
    input_values = list()
    result_ids = list()
    expected_results = list()

    for flag_str, value, expected_result, delay_ms in steps:
        flags = InferRequestHeader.FLAG_NONE
        if flag_str is not None:
            if "start" in flag_str:
                flags = flags | InferRequestHeader.FLAG_SEQUENCE_START
            if "end" in flag_str:
                flags = flags | InferRequestHeader.FLAG_SEQUENCE_END

        input_list = list()
        for b in range(batch_size):
            if input_dtype == np.object:
                in0 = np.full(tensor_shape, value, dtype=np.int32)
                in0n = np.array([str(x) for x in in0.reshape(in0.size)], dtype=object)
                in0 = in0n.reshape(tensor_shape)
            else:
                in0 = np.full(tensor_shape, value, dtype=input_dtype)
            input_list.append(in0)

        input_values.append(value)
        expected_results.append(expected_result)
        result_ids.append(ctx.async_run(
            { "INPUT" : input_list }, { "OUTPUT" : InferContext.ResultFormat.RAW},
            batch_size=batch_size, flags=flags))

        if delay_ms is not None:
            time.sleep(delay_ms / 1000.0)

    # Wait for the results in the order sent
    result = None
    for id, expected, value in zip(result_ids, expected_results, input_values):
        results = None
        while results == None:
            results = ctx.get_async_run_results(id, False)
            if results == None:
                if timeout_ms != None:
                    now_ms = int(round(time.time() * 1000))
                    if (now_ms - seq_start_ms) > timeout_ms:
                        raise TimeoutException("Timeout expired for {}".format(sequence_name))
                time.sleep(10.0 / 1000.0) # 10ms

        assert len(results) == 1
        assert "OUTPUT" in results
        result = results["OUTPUT"][0][0]
        if FLAGS.verbose:
            print("{} {}: + {} = {}".format(sequence_name, ctx.correlation_id(), value, result))

        if expected is not None:
            if input_dtype == np.object:
                assert int(result) == expected, "{}: expected result {}, got {}".format(
                    sequence_name, expected, int(result))
            else:
                assert result == expected, "{}: expected result {}, got {}".format(
                    sequence_name, expected, result)

def get_datatype(trial):
    # Get the datatype to use based on what models are available (see test.sh)
    if ("plan" in trial) or ("savedmodel" in trial):
        return np.float32
    if "graphdef" in trial:
        return np.dtype(object)
    return np.int32

def sequence_valid(ctx, rng, trial, model_name, dtype, len_mean, len_stddev, sequence_name):
    # Create a variable length sequence with "start" and "end" flags.
    seqlen = max(1, int(rng.normal(len_mean, len_stddev)))
    print("{} {}: valid seqlen = {}".format(sequence_name, ctx.correlation_id(), seqlen))

    values = rng.randint(0, 1024*1024, size=seqlen, dtype=dtype)

    steps = []
    expected_result = 0

    for idx, step in enumerate(range(seqlen)):
        flags = ""
        if idx == 0:
            flags += ",start"
        if idx == (seqlen - 1):
            flags += ",end"

        val = values[idx]
        delay_ms = None
        expected_result += val

        # (flag_str, value, expected_result, delay_ms)
        steps.append((flags, val, expected_result, delay_ms),)

    check_sequence_async(ctx, trial, model_name, dtype, steps,
                         sequence_name=sequence_name)

def sequence_valid_valid(ctx, rng, trial, model_name, dtype, len_mean, len_stddev, sequence_name):
    # Create two variable length sequences with "start" and "end"
    # flags, where both sequences use the same correlation ID and are
    # sent back-to-back.
    seqlen = [ max(1, int(rng.normal(len_mean, len_stddev))),
               max(1, int(rng.normal(len_mean, len_stddev))) ]
    print("{} {}: valid-valid seqlen[0] = {}, seqlen[1] = {}".format(
        sequence_name, ctx.correlation_id(), seqlen[0], seqlen[1]))

    values = [ rng.randint(0, 1024*1024, size=seqlen[0], dtype=dtype),
               rng.randint(0, 1024*1024, size=seqlen[1], dtype=dtype) ]

    for p in [0, 1]:
        steps = []
        expected_result = 0

        for idx, step in enumerate(range(seqlen[p])):
            flags = ""
            if idx == 0:
                flags += ",start"
            if idx == (seqlen[p] - 1):
                flags += ",end"

            val = values[p][idx]
            delay_ms = None
            expected_result += val

            # (flag_str, value, expected_result, delay_ms)
            steps.append((flags, val, expected_result, delay_ms),)

    check_sequence_async(ctx, trial, model_name, dtype, steps,
                         sequence_name=sequence_name)

def sequence_valid_no_end(ctx, rng, trial, model_name, dtype, len_mean, len_stddev, sequence_name):
    # Create two variable length sequences, the first with "start" and
    # "end" flags and the second with no "end" flag, where both
    # sequences use the same correlation ID and are sent back-to-back.
    seqlen = [ max(1, int(rng.normal(len_mean, len_stddev))),
               max(1, int(rng.normal(len_mean, len_stddev))) ]
    print("{} {}: valid-no-end seqlen[0] = {}, seqlen[1] = {}".format(
        sequence_name, ctx.correlation_id(), seqlen[0], seqlen[1]))

    values = [ rng.randint(0, 1024*1024, size=seqlen[0], dtype=dtype),
               rng.randint(0, 1024*1024, size=seqlen[1], dtype=dtype) ]

    for p in [0, 1]:
        steps = []
        expected_result = 0

        for idx, step in enumerate(range(seqlen[p])):
            flags = ""
            if idx == 0:
                flags += ",start"
            if (p == 0) and (idx == (seqlen[p] - 1)):
                flags += ",end"

            val = values[p][idx]
            delay_ms = None
            expected_result += val

            # (flag_str, value, expected_result, delay_ms)
            steps.append((flags, val, expected_result, delay_ms),)

    check_sequence_async(ctx, trial, model_name, dtype, steps,
                         sequence_name=sequence_name)

def sequence_no_start(ctx, rng, trial, model_name, dtype, sequence_name):
    # Create a sequence without a "start" flag. Sequence should get an
    # error from the server.
    seqlen = 1
    print("{} {}: no-start seqlen = {}".format(sequence_name, ctx.correlation_id(), seqlen))

    values = rng.randint(0, 1024*1024, size=seqlen, dtype=dtype)

    steps = []

    for idx, step in enumerate(range(seqlen)):
        flags = None
        val = values[idx]
        delay_ms = None

        # (flag_str, value, expected_result, delay_ms)
        steps.append((flags, val, None, delay_ms),)

    try:
        check_sequence_async(ctx, trial, model_name, dtype, steps,
                             sequence_name=sequence_name)
        assert False, "expected inference failure from missing START flag"
    except InferenceServerException as ex:
        if "must specify the START flag" not in ex.message():
            raise

def sequence_no_end(ctx, rng, trial, model_name, dtype, len_mean, len_stddev, sequence_name):
    # Create a variable length sequence with "start" flag but that
    # never ends. The sequence should be aborted by the server and its
    # slot reused for another sequence.
    seqlen = max(1, int(rng.normal(len_mean, len_stddev)))
    print("{} {}: no-end seqlen = {}".format(sequence_name, ctx.correlation_id(), seqlen))

    values = rng.randint(0, 1024*1024, size=seqlen, dtype=dtype)

    steps = []
    expected_result = 0

    for idx, step in enumerate(range(seqlen)):
        flags = ""
        if idx == 0:
            flags = "start"

        val = values[idx]
        delay_ms = None
        expected_result += val

        # (flag_str, value, expected_result, delay_ms)
        steps.append((flags, val, expected_result, delay_ms),)

    check_sequence_async(ctx, trial, model_name, dtype, steps,
                         sequence_name=sequence_name)

def stress_thread(name, seed, pass_cnt, correlation_id_base, trial, model_name, dtype):
    # Thread responsible for generating sequences of inference
    # requests.
    global _thread_exceptions

    print("Starting thread {} with seed {}".format(name, seed))
    rng = np.random.RandomState(seed)

    try:
        # Must use streaming GRPC context to ensure each sequences'
        # requests are received in order. Create 2 common-use contexts
        # with different correlation IDs that are used for most
        # inference requests. Also create some rare-use contexts that
        # are used to make requests with rarely-used correlation IDs.
        #
        # Need to remember the last choice for each context since we
        # don't want some choices to follow others since that gives
        # results not expected. See below for details.
        common_cnt = 2
        rare_cnt = 8
        ctxs = []
        last_choices = []

        for c in range(common_cnt + rare_cnt):
            ctxs.append(
                InferContext("localhost:8001", ProtocolType.GRPC, model_name,
                             correlation_id=correlation_id_base + c, streaming=True,
                             verbose=FLAGS.verbose))
            last_choices.append(None)

        rare_idx = 0
        for p in range(pass_cnt):
            # Common or rare context?
            if rng.rand() < 0.1:
                # Rare context...
                choice = rng.rand()
                ctx_idx = common_cnt + rare_idx

                # Send a no-end, valid-no-end or valid-valid
                # sequence... because it is a rare context this should
                # exercise the idle sequence path of the sequence
                # scheduler
                if choice < 0.33:
                    sequence_no_end(ctxs[ctx_idx], rng, trial, model_name, dtype,
                                    SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                    sequence_name=name)
                    last_choices[ctx_idx] = "no-end"
                elif choice < 0.66:
                    sequence_valid_no_end(ctxs[ctx_idx], rng, trial, model_name, dtype,
                                   SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                   sequence_name=name)
                    last_choices[ctx_idx] = "valid-no-end"
                else:
                    sequence_valid_valid(ctxs[ctx_idx], rng, trial, model_name, dtype,
                                   SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                   sequence_name=name)
                    last_choices[ctx_idx] = "valid-valid"

                rare_idx = (rare_idx + 1) % rare_cnt
            else:
                # Common context...
                ctx_idx = 0 if rng.rand() < 0.5 else 1
                ctx = ctxs[ctx_idx]
                last_choice = last_choices[ctx_idx]

                choice = rng.rand()

                # no-start cannot follow no-end since the server will
                # just assume that the no-start is a continuation of
                # the no-end sequence instead of being a sequence
                # missing start flag.
                if ((last_choice != "no-end") and
                    (last_choice != "valid-no-end") and
                    (choice < 0.01)):
                    sequence_no_start(ctx, rng, trial, model_name, dtype,
                                      sequence_name=name)
                    last_choices[ctx_idx] = "no-start"
                elif choice < 0.05:
                    sequence_no_end(ctx, rng, trial, model_name, dtype,
                                    SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                    sequence_name=name)
                    last_choices[ctx_idx] = "no-end"
                elif choice < 0.10:
                    sequence_valid_no_end(ctx, rng, trial, model_name, dtype,
                                   SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                   sequence_name=name)
                    last_choices[ctx_idx] = "valid-no-end"
                elif choice < 0.15:
                    sequence_valid_valid(ctx, rng, trial, model_name, dtype,
                                   SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                   sequence_name=name)
                    last_choices[ctx_idx] = "valid-valid"
                else:
                    sequence_valid(ctx, rng, trial, model_name, dtype,
                                   SEQUENCE_LENGTH_MEAN, SEQUENCE_LENGTH_STDEV,
                                   sequence_name=name)
                    last_choices[ctx_idx] = "valid"

    except Exception as ex:
        _thread_exceptions_mutex.acquire()
        try:
            _thread_exceptions.append(traceback.format_exc())
        finally:
            _thread_exceptions_mutex.release()
    print("Exiting thread {}".format(name))

def check_status(model_name):
    ctx = ServerStatusContext("localhost:8000", ProtocolType.HTTP, model_name, FLAGS.verbose)
    ss = ctx.get_server_status()
    print(ss)

if __name__ == '__main__':
    parser = argparse.ArgumentParser()
    parser.add_argument('-v', '--verbose', action="store_true", required=False, default=False,
                        help='Enable verbose output')
    parser.add_argument('-r', '--random-seed', type=int, required=False,
                        help='Random seed.')
    parser.add_argument('-t', '--concurrency', type=int, required=False, default=8,
                        help='Request concurrency. Default is 8.')
    parser.add_argument('-i', '--iterations', type=int, required=False, default=200,
                        help='Number of iterations of stress test to run. Default is 200.')
    FLAGS = parser.parse_args()

    # Initialize the random seed. For reproducibility each thread
    # maintains its own RNG which is initialized based on this seed.
    randseed = 0
    if FLAGS.random_seed != None:
        randseed = FLAGS.random_seed
    else:
        randseed = int(time.time())
    np.random.seed(randseed)

    print("random seed = {}".format(randseed))
    print("concurrency = {}".format(FLAGS.concurrency))
    print("iterations = {}".format(FLAGS.iterations))

    trial = "custom"
    dtype = get_datatype(trial)
    model_name = tu.get_sequence_model_name(trial, dtype)

    threads = []
    for idx, thd in enumerate(range(FLAGS.concurrency)):
        thread_name = "thread_{}".format(idx)

        # Create the seed for the thread. Since these are created in
        # reproducible order off of the initial seed we will get
        # reproducible results when given the same seed.
        seed = np.random.randint(2**32)

        # Each thread is reserved a block of correlation IDs or size
        # CORRELATION_ID_BLOCK_SIZE
        correlation_id_base = 1 + (idx * CORRELATION_ID_BLOCK_SIZE)

        threads.append(threading.Thread(
                    target=stress_thread,
                    args=(thread_name, seed, FLAGS.iterations,
                          correlation_id_base, trial, model_name, dtype)))

    for t in threads:
        t.start()
    for t in threads:
        t.join()

    check_status(model_name)

    _thread_exceptions_mutex.acquire()
    try:
        if len(_thread_exceptions) > 0:
            for ex in _thread_exceptions:
                print("*********\n{}".format(ex))
            sys.exit(1)
    finally:
        _thread_exceptions_mutex.release()

    sys.exit(0)

max_stars_repo_name: AchimLoerke/tensorrt-inference-server

max_issues_repo_name: AchimLoerke/tensorrt-inference-server

max_forks_repo_name: AchimLoerke/tensorrt-inference-server


text:
import h5py
import numpy as np
import GMatElastoPlastic.Cartesian3d as GMat

with h5py.File('Cartesian3d_random.hdf5', 'w') as data:

    nelem = 1000
    nip = 4
    iden = 2.0 * np.random.random([nelem, nip])
    iden = np.where(iden < 1.0, 0.0, iden)
    iden = np.where(iden >= 1.0, 1.0, iden)
    iden = iden.astype(np.int)

    shape = np.array([nelem, nip], np.int)

    data['/shape'] = shape

    mat = GMat.Array2d(shape)

    I = np.where(iden == 0, 1, 0).astype(np.int)
    n = np.sum(I)
    idx = np.zeros(I.size, np.int)
    idx[np.argwhere(I.ravel() == 1).ravel()] = np.arange(n)
    idx = idx.reshape(I.shape)
    K = np.ones(n)
    G = np.ones(n)
    sigy0 = 0.1 * np.ones(n)
    H = np.ones(n)

    data['/LinearHardening/I'] = I
    data['/LinearHardening/idx'] = idx
    data['/LinearHardening/K'] = K
    data['/LinearHardening/G'] = G
    data['/LinearHardening/sigy0'] = sigy0
    data['/LinearHardening/H'] = H

    mat.setLinearHardening(I, idx, K, G, sigy0, H)

    I = np.where(iden == 1, 1, 0).astype(np.int)
    n = np.sum(I)
    idx = np.zeros(I.size, np.int)
    idx[np.argwhere(I.ravel() == 1).ravel()] = np.arange(n)
    idx = idx.reshape(I.shape)
    K = np.ones(n)
    G = np.ones(n)

    data['/Elastic/I'] = I
    data['/Elastic/idx'] = idx
    data['/Elastic/K'] = K
    data['/Elastic/G'] = G

    mat.setElastic(I, idx, K, G)

    for i in range(20):

        mat.increment()

        GradU = 20 * np.random.random([nelem, nip, 3, 3])

        data['/random/{0:d}/GradU'.format(i)] = GradU

        Eps = np.einsum('...ijkl,...lk->...ij', mat.I4s(), GradU)
        mat.setStrain(Eps)

        data['/random/{0:d}/Stress'.format(i)] = mat.Stress()
        data['/random/{0:d}/Tangent'.format(i)] = mat.Tangent()
        data['/random/{0:d}/Epsp'.format(i)] = mat.Epsp()

max_stars_repo_name: tdegeus/GMatElastoPlastic

max_issues_repo_name: tdegeus/GMatElastoPlastic

max_forks_repo_name: tdegeus/GMatElastoPlastic


text:
#!/usr/bin/env python

from __future__ import print_function

import numpy as np
import sys
import benchmark_decorator as dectimer

#--------------------
# Function: integrand
#--------------------
integrand = lambda x: np.exp(x)

#-----------------------------
# Function: compute_quadrature
#-----------------------------
@dectimer.bench_time(3)
def compute_quadrature(n):
    """
      Perform the Gauss-Legendre Quadrature at the prescribed order n
    """
    a = -3.0
    b = 3.0

    # Gauss-Legendre (default interval is [-1, 1])
    x, w = np.polynomial.legendre.leggauss(n)

    # Translate x values from the interval [-1, 1] to [a, b]
    t = 0.5*(x + 1)*(b - a) + a

    gauss = sum(w * integrand(t)) * 0.5*(b - a)

if len(sys.argv) < 1:
    print('Usage:')
    print('     python ' + sys.argv[0] + ' N')
    print('Please specify the order of the quadrature.')
    sys.exit()

order = int(sys.argv[1])

print('Gauss-Legendre Quadrature of order: ', order)

compute_quadrature(order)

print(' ')


text:
# Function for testing models

import numpy as np
import torch

def test(data, model, pad_idx, device, args):
    model.eval()
    with torch.no_grad():
        all_correct_trials = [] # list of booleans indicating whether correct
        for batch in data:
            out, attn_wts = model(batch.src, batch.trg)
            preds = torch.argmax(out,dim=2)
            correct_pred = preds == batch.trg
            correct_pred = correct_pred.cpu().numpy()
            mask = batch.trg == pad_idx # mask out padding
            mask = mask.cpu().numpy()
            correct = np.logical_or(mask,correct_pred)
            correct = correct.all(0).tolist()
            all_correct_trials += correct

    accuracy = np.mean(all_correct_trials)
    model.train()
    return accuracy

max_issues_repo_name: JanAthmer/Scriptie

max_forks_repo_name: JanAthmer/Scriptie


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-
"""
Created on Mon May  6 13:50:57 2019

@author: jishnu
"""
#%%
import pandas as pd
import numpy as np
import matplotlib.pyplot as plt
from astropy import time, coordinates as coord, units as u
from astropy.io import fits
from astropy.time import Time
import glob
import os

#%%
# =============================================================================
# This script runs through files containing detected sources from Sextractor and Psfex
# source files have the following columns, refer back to photometry script
# and sextractor to add/remove other columns
# Index(['X_IMAGE', 'Y_IMAGE', 'ALPHAWIN_J2000', 'DELTAWIN_J2000', 'SNR_WIN',
#        'MAG_POINTSOURCE', 'MAGERR_POINTSOURCE'],
#       dtype='object')
#
# The output is a dataframe and Lightcurve csv written to disk
# =============================================================================
#%%
'''
The coordinates for variable and check stars
'''
V_x,V_y = 1230,1750
C1x,C2x = 1670,288
C2x,C2y = 960,2054

#%%
#root_folder = '../../'
#os.chdir(root_folder)
os.chdir('/home/jishnu/JCBT/TYC3315/14Nov2/reduced/sources')
os.chdir('..')
ffiles = glob.glob('*proc.fits')
ffiles.sort()

#%%
def read(filename): #read a csv file and returns a dataframe
    df = pd.read_csv(filename,sep='\t')
    df = df.drop([df.columns[0]],axis=1)
    return df

#%%
def search(x,y,df,t):
    '''
    #here x,y are pixel coords of target star, t is the size of the box
    used for search
    '''
    flt = df[ #setting condtion so to consider sources inside a box around(x,y)
            (df['X_IMAGE'] < (x+t)) &
            (df['X_IMAGE'] > (x-t)) &
            (df['Y_IMAGE'] < (x+t)) &
            (df['Y_IMAGE'] > (x-t))
            ]
    print(flt)
    return flt['MAG_POINTSOURCE'],flt['MAGERR_POINTSOURCE'] #returns pandas series object

#%%
def cone_search(x,y,df,dist):
    '''
    #another implematation of previous function
    takes in (x,y) coordinates, a dataframe and search radius
    returns the magnitude and error if any source is
    found around the given coordinates
    '''
    #uses distance formula to check how close the point is, to given threshold
    fdf = df[
            (np.sqrt((df['X_IMAGE']-x)**2 +(df['Y_IMAGE']-y)**2) < dist)
            ]
    print(fdf)

    return fdf['MAG_POINTSOURCE'],fdf['MAGERR_POINTSOURCE']

#%%
def coord_search(ra,dec,df,dist):
    '''
    Uses RA,DEC instead of pixel coordinates to search
    Use this function only if there is accurate WCS information and the
    fits file has been platesolved
    '''
    fdf = df[ #uses distance formula to check how close the point is, to given threshold
            (np.sqrt((df['ALPHAWIN_J2000']-ra)**2 +
                     (df['DELTAWIN_J2000']-dec)**2) < dist)
            ]
    print(fdf)
    if len(fdf['MAG_POINTSOURCE'])==0:
        return 0,0 #returns 0 if no sources found
    else:
        return fdf['MAG_POINTSOURCE'],fdf['MAGERR_POINTSOURCE']
#%%
def sorting(x,y,files,t):
    '''
    Takes (x,y), files and distance threshold and returns the magnitude
    and error of the source located at (x,y)
    '''
    mag,err = [],[]
    for f in files: #honestly this could have been done in one function..
        print(f+'\n'*2)
        mags = cone_search(x,y,read(f),t) #pick your poison here
        mag.append(mags[0])
        err.append(mags[1])
        print('\n'*2+'-'*50)
    return mag,err #returns the magnitude and error columns, as pandas series

#%%
def filter(lst,err): #serves same purpose as the above function, now you might ask why there are two tools for the same job, to which i have no answer. They are nice, first one i came up with, in my head; and second one, well Internet people! https://stackoverflow.com/questions/12141150/from-list-of-integers-get-number-closest-to-a-given-value
    '''
    Input is pandas series of mag and error, if there are more than one elements
    go into the function, else just return the value from series
    Takes the magnitude and associated error columns, and removes duplicates
    (sometimes there are more than 1 source in the given (x,y) and d search region)
    The function tries to find the correct source by taking mean of the list
    and looking to reject other sources based on that

    Should not be used if the sources are clumped up and has the same mag

    Needs better filtering, until then; stuck with this dirty way
    '''
    opt,Err = [0],[]
    for i in range(len(lst)):
        try:
            if len(lst[i])==1:
                opt.append(lst[i].iloc[0])
                Err.append(err[i].iloc[0])
            else: #fancy pants way to do the earlier function
                opt.append(min([lst[i].iloc[n] for n in range(len(lst[i]))],
                               key=lambda x:abs(x-np.mean(opt))))
                Err.append(min([err[i].iloc[n] for n in range(len(err[i]))],
                               key=lambda x:abs(x-np.mean(opt))))
        except:
            opt.append(0)
            Err.append(0)

    return opt[1:],Err #returns 2 lists, mag and err. Do what you want with them

#%%
def Filter(lst,err):
    opt = [0] # creating list with one element to avoid np raising error
              # while taking mean
    Err = []
    for i in range(len(lst)):
        if len(lst[i])==1:
            opt.append(lst[i].iloc[0])
            Err.append(err[i].iloc[0])
        else:
            if (abs(np.mean(opt)-lst[i].iloc[0])) > (abs(np.mean(opt)-lst[i].iloc[1])):
                opt.append(lst[i].iloc[1])
                Err.append(err[i].iloc[1])
            else:
                opt.append(lst[i].iloc[0])
                Err.append(err[i].iloc[0])

    return opt[1:],Err

#%%
def get_info(imageName):
    loc = coord.EarthLocation.of_site('vbo')
    hdul = fits.open(imageName)
    ut=(hdul[0].header['UT'])
    ra=(hdul[0].header['RA'])
    dec=(hdul[0].header['DEC'])
    am=(hdul[0].header['AIRMASS'])
    t = Time(hdul[0].header['DATE-OBS'], scale='utc')
    jd=(t.jd) #we only use JD.

    obs_target = coord.SkyCoord(hdul[0].header['RA'],hdul[0].header['DEC'],
                                unit = (u.hourangle,u.deg),frame = 'icrs'
                                )
    times = time.Time(t.jd, format='mjd',scale='utc', location=loc
                      )
    ltt_bary = times.light_travel_time(obs_target)
    time_barycentre = times.tdb + ltt_bary
    bjd=(time_barycentre.value)

    hdul.close()
    return ut,ra,dec,am,jd,bjd
    #USAGE
    #ut,ra,dec,am,jd,bjd = get_info(imagename)
#%%
df = pd.DataFrame()
ut,ra,dec,am,jd,bjd = [],[],[],[],[],[]
#%%
def main():
    for f in range(len(ffiles)):
        UT,RA,DEC,AM,JD,BJD = get_info(ffiles[f])
        ut.append(UT)
        ra.append(RA)
        dec.append(DEC)
        am.append(AM)
        jd.append(JD)
        bjd.append(BJD)

    folder = 'sources/'
    os.chdir(folder)

    Files = [i+'_source_cat.csv' for i in ffiles] #can just read the files off disk but where's the fun in that.
    m,e = sorting(V_x,V_y,Files,100) #loopy logic of this function! yay me
    Vmag,Verr = filter(m,e)

    m,e = sorting(C1x,C2x,Files,100)
    C1mag,C1err = filter(m,e)

    m,e = sorting(C2x,C2y,Files,100)
    C2mag,C2err = filter(m,e)

#    m,e = sorting(50.381181,47.433372,Files,0.01) #loopy logic of this function! yay me
#    Vmag,Verr = filter(m,e)

    df = pd.DataFrame()
    df = df.assign(
            **{'UT':ut,'RA':ra,'DEC':dec,'airmass':am,
               'JD':jd,'BJD':bjd,'Vmag':Vmag,'v_err':Verr,
               'C1':C1mag,'C1_err':C1err,'C2':C2mag,'C2_err':C2err
                    })
    df['Vmag_d'] = df['Vmag']-df['C2']
    df.to_csv('LC.csv')
    return df

#%%
plt.figure(figsize=(8,6))
plt.plot(bjd,df.Vmag,'o',label='v')
plt.plot(bjd,df.C1,'o',label='C1')
plt.plot(bjd,df.C2,'o',label='C2')
plt.legend()

#%%
plt.plot(bjd,df.Vmag-df.C1,'o',label='v-c1')
plt.plot(bjd,df.Vmag-df.C2,'o',label='v-c2')
plt.plot(bjd,df.C1-df.C2,'o',label='c1-c2')
plt.legend()

max_stars_repo_name: jishnupdas/Photometry-tools

max_issues_repo_name: jishnupdas/Photometry-tools

max_forks_repo_name: jishnupdas/Photometry-tools

reason: import numpy,from astropy


text:
# %%
import os

import matplotlib.pyplot as plt
import numpy as np
import pandas as pd
import torch
from ml_matrics import density_scatter, precision_recall_curve, roc_curve
from sklearn.model_selection import KFold
from torch.utils.data import DataLoader
from torch.utils.tensorboard import SummaryWriter

from atom_joggling.cgcnn import CGCNN, CrystalGraphData, collate_batch
from atom_joggling.utils import ROOT, mean, parser

# %%
args, _ = parser.parse_known_args()

data_df = pd.read_csv(f"{ROOT}/{args.data_path}").iloc[:, 0:3]

n_splits = 5
kfold = KFold(n_splits, random_state=0, shuffle=True)

untrained_perfs, results = [], []

print(f"\n- Task: {(task := args.task)}")
print(f"- Using CUDA: {(use_cuda := torch.cuda.is_available())}")
print(f"- Batch size: {(batch_size := args.batch_size):,d}")
print(f"- Epochs: {(epochs := args.epochs):,d}")
print(f"- Number of samples: {len(data_df):,d}")
print(f"- Verbose: {(verbose := args.verbose)}")
print(f"- Joggle: {(joggle := args.joggle)}")

csv_name = os.path.splitext(args.data_path)[0].split(os.sep)[-1]
default_dir = f"{ROOT}/runs/cgcnn/{task}"
default_dir += f"/{csv_name}-{n_splits}folds-{epochs}epochs-{batch_size}batch"
print(f"- Output directory: {(out_dir := args.out_dir or default_dir)}\n")

loader_args = {
    "batch_size": batch_size,
    "collate_fn": lambda batch: collate_batch(batch, use_cuda),
}

targets = data_df.iloc[:, 2]
if task == "classification":
    print(f"dummy accuracy: {targets.mean():.3f}")
else:  # regression
    print(f"dummy MAE: {(targets - targets.mean()).abs().mean():.3f}")

# %%
for fold, (train_idx, test_idx) in enumerate(kfold.split(data_df), 1):
    print(f"\nFold {fold}/{n_splits}")
    fold_dir = f"{out_dir}/fold_{fold}"
    train_df, test_df = data_df.iloc[train_idx], data_df.iloc[test_idx]

    train_set = CrystalGraphData(train_df, task, joggle=joggle)
    test_set = CrystalGraphData(test_df, task)

    model = CGCNN(
        task,
        args.robust,
        train_set.elem_emb_len,
        train_set.nbr_fea_len,
        n_targets=train_set.n_targets,
        checkpoint_dir=fold_dir,
    )
    if use_cuda:
        model.cuda()

    optimizer = torch.optim.AdamW(model.parameters())

    train_loader = DataLoader(train_set, **loader_args)
    test_loader = DataLoader(test_set, **loader_args)

    writer = SummaryWriter(fold_dir)

    # sanity check: test untrained model performance
    material_ids, formulas, targets, outputs = model.predict(test_loader)

    if task == "classification":
        untrained_acc = (targets == outputs.argmax(dim=1)).float().mean()
        untrained_perfs.append(untrained_acc)
        print(f"untrained accuracy: {untrained_acc:.3f}")
    else:  # regression
        untrained_mae = (targets - outputs).abs().mean()
        untrained_perfs.append(untrained_mae)
        print(f"untrained MAE: {untrained_mae:.3f}")

    model.fit(train_loader, test_loader, optimizer, epochs, writer, verbose=verbose)

    mp_ids, formulas, targets, preds = model.predict(test_loader)

    out_df = pd.DataFrame(  # collect model output into dataframe
        zip([fold] * len(mp_ids), mp_ids, formulas, targets.numpy()),
        columns=["fold", "material_id", "formula", "target"],
    )

    if task == "classification":
        out_df["pred"] = preds.argmax(dim=1)
        for idx, arr in enumerate(preds.softmax(dim=1).numpy().T):
            out_df[f"softmax_{idx}"] = arr

        test_acc = (targets == preds.argmax(dim=1)).float().mean()
        results.append([out_df, test_acc])
        print(f"test accuracy: {test_acc:.3f}")

    else:  # regression
        if model.robust:
            preds, log_std = preds.chunk(2, dim=1)
            out_df["std"] = log_std.exp()
        out_df["pred"] = preds

        test_mae = (targets - preds).abs().mean()
        results.append([out_df, test_mae])
        print(f"test MAE: {test_mae:.3f}")

# %%
untrained_perf_str = f"{mean(untrained_perfs):.3g} +/- {np.std(untrained_perfs):.3g}"
out_dfs, test_accs = zip(*results)

out_df = pd.concat(out_dfs)
out_df.to_csv(f"{out_dir}/output.csv", index=False)

test_perf_str = f"{mean(test_accs):.3g} +/- {np.std(test_accs):.3g}"

print(f"\n\n{n_splits}-fold split results")

if task == "classification":

    if len(out_df) < 100:
        out_df.plot.bar(x="formula", y=["softmax_1", "target"], figsize=[18, 8])
        plt.savefig(out_dir + "/cgcnn_val_preds.png", dpi=200, bbox_inches="tight")
        plt.close()

    roc_auc, _ = roc_curve(out_df.target, out_df.softmax_1)
    plt.savefig(out_dir + "/roc_auc_curve.png", dpi=200, bbox_inches="tight")
    plt.close()
    prec, _ = precision_recall_curve(out_df.target, out_df.softmax_1)
    plt.savefig(out_dir + "/precision_recall_curve.png", dpi=200, bbox_inches="tight")
    plt.close()

    print(f"mean untrained accuracy: {untrained_perf_str}")
    print(f"mean test accuracy: {test_perf_str}")
    print(f"mean ROC AUC: {roc_auc:.3g}")
    print(f"mean precisions: {prec:.3g}")

else:  # regression
    print(f"mean untrained MAE: {untrained_perf_str}")
    print(f"mean test MAE: {test_perf_str}")

    density_scatter(out_df.target, out_df.pred)

max_stars_repo_name: janosh/atom-joggling

max_issues_repo_name: janosh/atom-joggling

max_forks_repo_name: janosh/atom-joggling


text:
import numpy as np
from scipy.optimize import minimize

from optimize_utils import *

def _minimize_rhc(*args, **kwargs):
    # randomized hill-climbing
    options = kwargs["options"]
    method = options.pop("method")
    kwargs["method"] = method
    remaining = options["maxiter"]
    best_result = None
    while remaining > 0:
        options["maxiter"] = remaining
        result = minimize(*args, **kwargs)
        if best_result is None or best_result.fun > result.fun:
            best_result = result
        remaining -= result.nfev
    return best_result

def _minimize_twiddle(fun, x0, args=(), callback=None, tol=0.2, maxiter=1000):
    nfev = 0
    def score(x):
        # nfev += 1  # TODO uncomment when nonlocal keyword can be used
        res = fun(x, *args)
        if callback is not None:
            callback(x)
        return res

    p = x0.copy()
    num_params = len(p)
    dp = np.ones(num_params)
    best_err = score(p)
    nfev += 1
    stop = False
    n_iter = 0
    while dp.sum() > tol and not stop:
        n_iter += 1
        stop = True
        for i in range(num_params):
            p[i] += dp[i]
            err = score(p)
            nfev += 1
            if err < best_err:
                best_err = err
                dp[i] *= 1.1
            else:
                p[i] -= 2 * dp[i]
                err = score(p)
                nfev += 1
                if err < best_err:
                    best_err = err
                    dp[i] *= 1.1
                else:
                    p[i] += dp[i]
                    dp[i] *= 0.9
            # subtraction 1 since each loop is 2 iterations
            if nfev >= maxiter - 1:
                break
        else:
            stop = False
    return to_result(x=p, fun=best_err, niter=n_iter, nfev=nfev)