text:
# You may obtain a copy of the License at
#
#     http://www.apache.org/licenses/LICENSE-2.0
#
# Unless required by applicable law or agreed to in writing, software
# distributed under the License is distributed on an "AS IS" BASIS,
# WITHOUT WARRANTIES OR CONDITIONS OF ANY KIND, either express or implied.
# See the License for the specific language governing permissions and
# limitations under the License.
# ==============================================================================
"""Utilities for building transformer layers."""
import math
import numpy as np
import tensorflow as tf
from tensorflow.python.util import nest

def log_prob_from_logits(logits, reduce_axis=-1):
  """return log prob use log sum func"""
  return logits - tf.reduce_logsumexp(logits, axis=reduce_axis, keepdims=True)

#
def shape_list(tensor):
  """Return list of dims, statically where possible."""
  tensor = tf.convert_to_tensor(tensor)

  if tensor.get_shape().dims is None:
    return tf.shape(tensor)

  static = tensor.get_shape().as_list()
  shape = tf.shape(tensor)

  ret = []
  for i, _ in enumerate(static):
    dim = static[i]
    if dim is None:
      dim = shape[i]
    ret.append(dim)
  return ret

def merge_beam_dim(tensor):
  """Reshapes first two dimensions in to single dimension."""
  shape = shape_list(tensor)
  shape[0] *= shape[1]  # batch -> batch * beam_size
  shape.pop(1)  # Remove beam dim
  return tf.reshape(tensor, shape)

def unmerge_beam_dim(tensor, batch_size, beam_size):
  """Reshapes first dimension back to [batch_size, beam_size]."""
  shape = shape_list(tensor)
  new_shape = [batch_size] + [beam_size] + shape[1:]
  return tf.reshape(tensor, new_shape)

def expand_to_beam_size(tensor, beam_size):
  """Tiles a given tensor by beam_size."""
  tensor = tf.expand_dims(tensor, axis=1)
  tile_dims = [1] * tensor.shape.ndims
  tile_dims[1] = beam_size

  return tf.tile(tensor, tile_dims)

def get_state_shape_invariants(tensor):
  """Returns the shape of the tensor but sets middle dims to None."""
  shape = tensor.shape.as_list()
  for i in range(1, len(shape) - 1):
    shape[i] = None
  return tf.TensorShape(shape)

def compute_batch_indices(batch_size, beam_size):
  """Computes the i'th coordinate that contains the batch index for gathers."""
  batch_pos = tf.range(batch_size * beam_size) // beam_size
  batch_pos = tf.reshape(batch_pos, [batch_size, beam_size])
  return batch_pos

def create_make_unique(inputs):
  """Replaces the lower bits of each element with iota."""
  if inputs.shape.ndims != 2:
    raise ValueError("Input of top_k_with_unique must be rank-2 "
                     "but got: %s" % inputs.shape)

  height = inputs.shape[0]
  width = inputs.shape[1]
  zeros = tf.zeros([height, width], dtype=tf.int32)

  log2_ceiling = int(math.ceil(math.log(int(width), 2)))
  next_power_of_two = 1 << log2_ceiling
  count_mask = ~(next_power_of_two - 1)
  count_mask_r0 = tf.constant(count_mask)
  count_mask_r2 = tf.fill([height, width], count_mask_r0)

  smallest_normal = 1 << 23
  smallest_normal_r0 = tf.constant(smallest_normal, dtype=tf.int32)
  smallest_normal_r2 = tf.fill([height, width], smallest_normal_r0)

  low_bit_mask = ~(1 << 31)
  low_bit_mask_r0 = tf.constant(low_bit_mask, dtype=tf.int32)
  low_bit_mask_r2 = tf.fill([height, width], low_bit_mask_r0)

  iota = tf.tile(
      tf.expand_dims(tf.range(width, dtype=tf.int32), 0), [height, 1])

  input_r2 = tf.bitcast(inputs, tf.int32)
  abs_r2 = tf.bitwise.bitwise_and(input_r2, low_bit_mask_r2)
  if_zero_r2 = tf.equal(abs_r2, zeros)
  smallest_normal_preserving_sign_r2 = tf.bitwise.bitwise_or(
      input_r2, smallest_normal_r2)
  input_no_zeros_r2 = tf.where(if_zero_r2, smallest_normal_preserving_sign_r2,
                               input_r2)

  and_r2 = tf.bitwise.bitwise_and(input_no_zeros_r2, count_mask_r2)
  or_r2 = tf.bitwise.bitwise_or(and_r2, iota)
  return tf.bitcast(or_r2, tf.float32)

def create_topk_unique(inputs, k):
  """Creates the top k values in sorted order with indices."""
  height = inputs.shape[0]
  width = inputs.shape[1]
  neg_inf_r0 = tf.constant(-np.inf, dtype=tf.float32)
  ones = tf.ones([height, width], dtype=tf.float32)
  neg_inf_r2 = ones * neg_inf_r0
  inputs = tf.where(tf.is_nan(inputs), neg_inf_r2, inputs)

  tmp = inputs
  topk_r2 = tf.zeros([height, k], dtype=tf.float32)
  for i in range(k):
    kth_order_statistic = tf.reduce_max(tmp, axis=1, keepdims=True)
    k_mask = tf.tile(
        tf.expand_dims(tf.equal(tf.range(k), tf.fill([k], i)), 0), [height, 1])
    topk_r2 = tf.where(k_mask, tf.tile(kth_order_statistic, [1, k]), topk_r2)
    ge_r2 = tf.greater_equal(inputs, tf.tile(kth_order_statistic, [1, width]))
    tmp = tf.where(ge_r2, neg_inf_r2, inputs)

  log2_ceiling = int(math.ceil(math.log(float(int(width)), 2)))
  next_power_of_two = 1 << log2_ceiling
  count_mask = next_power_of_two - 1
  mask_r0 = tf.constant(count_mask)
  mask_r2 = tf.fill([height, k], mask_r0)
  topk_r2_s32 = tf.bitcast(topk_r2, tf.int32)
  topk_indices_r2 = tf.bitwise.bitwise_and(topk_r2_s32, mask_r2)
  return topk_r2, topk_indices_r2

def top_k_with_unique(inputs, k):
  """Finds the values and indices of the k largests entries."""
  unique_inputs = create_make_unique(tf.cast(inputs, tf.float32))
  top_values, indices = create_topk_unique(unique_inputs, k)
  top_values = tf.cast(top_values, inputs.dtype)
  return top_values, indices

def compute_topk_scores_and_seq(sequences,
                                scores,
                                scores_to_gather,
                                flags,
                                beam_size,
                                batch_size,
                                prefix="default",
                                states_to_gather=None):
  """Given sequences and scores, will gather the top k=beam size sequences."""
  _, topk_indexes = tf.nn.top_k(scores, k=beam_size)
  batch_pos = compute_batch_indices(batch_size, beam_size)
  top_coordinates = tf.stack([batch_pos, topk_indexes], axis=2)

  def gather(tensor, name):
    return tf.gather_nd(tensor, top_coordinates, name=(prefix + name))

  topk_seq = gather(sequences, "_topk_seq")
  topk_flags = gather(flags, "_topk_flags")
  topk_gathered_scores = gather(scores_to_gather, "_topk_scores")
  if states_to_gather:
    topk_gathered_states = nest.map_structure(
        lambda state: gather(state, "_topk_states"), states_to_gather)
  else:
    topk_gathered_states = states_to_gather

  return topk_seq, topk_gathered_scores, topk_flags, topk_gathered_states

def create_padding_mask(seq):
  """
  create padding mask
  """
  seq = tf.cast(tf.math.equal(seq, 0), tf.float32)
  return seq[:, tf.newaxis, tf.newaxis, :]  # (batch_size, 1, 1, seq_len)
  # seq = tf.expand_dims(seq, 1)
  # seq = tf.expand_dims(seq, 1)
  # return seq

def create_masks(inp, tar):
  """
  create encode and decode mask
  """
  # encoder self-attention mask
  enc_padding_mask = create_padding_mask(inp)
  # decoder self-attention mask
  dec_padding_mask = create_padding_mask(inp)

  # decoder context mask
  look_ahead_mask = create_look_ahead_mask(tf.shape(tar)[1])
  dec_target_padding_mask = create_padding_mask(tar)
  combined_mask = tf.maximum(dec_target_padding_mask, look_ahead_mask)

  return enc_padding_mask, combined_mask, dec_padding_mask

def create_look_ahead_mask(tar):
  """
  create look ahead mask for decode mask
  """
  size = tf.shape(tar)[1]
  look_ahead_mask = 1 - tf.linalg.band_part(tf.ones((size, size)), -1, 0)
  dec_target_padding_mask = create_padding_mask(tar)
  combined_mask = tf.maximum(dec_target_padding_mask, look_ahead_mask)
  return combined_mask  # (batch_size, 1, seq_len, seq_len)

max_stars_repo_name: didichuxing/delta

max_issues_repo_name: didichuxing/delta

max_forks_repo_name: didichuxing/delta


text:
from ..backend import Backend
from ..result import UnfoldingResult
import ROOT
import numpy as np
import root_numpy
import math

class SVDBackend(Backend):
    def __init__(self, regularization_param, bins_min, bins_max):
        """
        Parameters
        ----------
        regularization_param : number
            Regularization parameter to use in the unfolding method
        bins_min : number
            Lower limit of the histograms' bins
        bins_max : number
            Upper limit of the histograms' bins
        """
        self.kreg = regularization_param
        self.bins_min = bins_min
        self.bins_max = bins_max

    def solve(self, data, statcov, xini, bini, R):
        bins_min = self.bins_min
        bins_max = self.bins_max
        # Transform the reponse matrix from probabilities to events
        # R = np.multiply(xini, R, where=xini != 0)
        # Transform python arrays to ROOT TH1D
        datar = ROOT.TH1D("data", "data", data.shape[0], bins_min, bins_max)
        root_numpy.array2hist(data, datar)
        xinir = ROOT.TH1D("xini", "xini", xini.shape[0], bins_min, bins_max)
        root_numpy.array2hist(xini, xinir)
        binir = ROOT.TH1D("bini", "bini", bini.shape[0], bins_min, bins_max)
        root_numpy.array2hist(bini, binir)
        Adet = ROOT.TH2D(
            "R",
            "R",
            xini.shape[0],
            bins_min,
            bins_max,
            bini.shape[0],
            bins_min,
            bins_max,
        )
        root_numpy.array2hist(R, Adet)
        # Compute covariance matrix assuming the amount of events in each bin
        # come from a Poisson distribution
        statcovr = ROOT.TH2D(
            "statcov",
            "statcov",
            data.shape[0],
            bins_min,
            bins_max,
            data.shape[0],
            bins_min,
            bins_max,
        )
        root_numpy.array2hist(statcov, statcovr)
        # Create TSVDUnfold object and initialise
        tsvdunf = ROOT.TSVDUnfold(datar, statcovr, binir, xinir, Adet)
        # It is possible to normalise unfolded spectrum to unit area
        tsvdunf.SetNormalize(ROOT.kFALSE)
        # Perform the unfolding with regularisation parameter kreg = self.kreg
        # - the larger kreg, the finer grained the unfolding, but the more fluctuations occur
        # - the smaller kreg, the stronger is the regularisation and the bias
        unfres = tsvdunf.Unfold(self.kreg)
        # Get the distribution of the d to cross check the regularization
        # - choose kreg to be the point where |d_i| stop being statistically significantly >>1
        ddist = tsvdunf.GetD()
        self.kreg_distribution = root_numpy.hist2array(ddist)
        # Get the distribution of the singular values
        svdist = tsvdunf.GetSV()
        # Compute the error matrix for the unfolded spectrum using toy MC
        # using the measured covariance matrix as input to generate the toys
        # 100 toys should usually be enough
        # The same method can be used for different covariance matrices separately.
        ustatcov = tsvdunf.GetUnfoldCovMatrix(statcovr, 100)
        # Now compute the error matrix on the unfolded distribution originating
        # from the finite detector matrix statistics
        uadetcov = tsvdunf.GetAdetCovMatrix(100)
        # Sum up the two (they are uncorrelated)
        ustatcov.Add(uadetcov)
        # Get the computed regularized covariance matrix (always corresponding to total uncertainty passed in constructor) and add uncertainties from finite MC statistics.
        utaucov = tsvdunf.GetXtau()
        utaucov.Add(uadetcov)
        # Get the computed inverse of the covariance matrix
        uinvcov = tsvdunf.GetXinv()
        # Convert result to a python array
        result = root_numpy.hist2array(unfres)
        error = root_numpy.hist2array(utaucov)
        # the std to plot is just:
        # np.sqrt(np.diagonal(error))
        return UnfoldingResult(result, error)

ext: py

reason: import numpy


text:
# License: Apache-2.0
import databricks.koalas as ks
import numpy as np
import pandas as pd
import pytest
from pandas.testing import assert_frame_equal
from pyspark.ml.classification import RandomForestClassifier as RFCSpark
from xgboost import XGBClassifier, XGBRFClassifier

from gators.feature_selection.select_from_models import SelectFromModels

ks.set_option("compute.default_index_type", "distributed-sequence")

@pytest.fixture
def data():
    X = pd.DataFrame(
        {
            "A": [22.0, 38.0, 26.0, 35.0, 35.0, 28.11, 54.0, 2.0, 27.0, 14.0],
            "B": [7.25, 71.28, 7.92, 53.1, 8.05, 8.46, 51.86, 21.08, 11.13, 30.07],
            "C": [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
        }
    )
    y = pd.Series([0, 1, 1, 1, 0, 0, 0, 0, 1, 1], name="TARGET")
    X_expected = X[["A", "B"]].copy()
    model1 = XGBClassifier(
        random_state=0,
        subsample=1.0,
        n_estimators=2,
        max_depth=2,
        eval_metric="logloss",
        use_label_encoder=False,
    )
    model2 = XGBRFClassifier(
        random_state=0,
        subsample=1.0,
        n_estimators=2,
        max_depth=2,
        eval_metric="logloss",
        use_label_encoder=False,
    )
    obj = SelectFromModels(models=[model1, model2], k=2).fit(X, y)
    return obj, X, X_expected

@pytest.fixture
def data_ks():
    X = ks.DataFrame(
        {
            "A": [22.0, 38.0, 26.0, 35.0, 35.0, 28.11, 54.0, 2.0, 27.0, 14.0],
            "B": [7.25, 71.28, 7.92, 53.1, 8.05, 8.46, 51.86, 21.08, 11.13, 30.07],
            "C": [0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0, 0.0],
        }
    )
    y = ks.Series([0, 1, 1, 1, 0, 0, 0, 0, 1, 1], name="TARGET")
    X_expected = X[["A", "B"]].to_pandas().copy()
    model1_ks = RFCSpark(numTrees=2, maxDepth=2, labelCol=y.name, seed=0)
    mode2_ks = RFCSpark(numTrees=1, maxDepth=2, labelCol=y.name, seed=0)
    obj = SelectFromModels(models=[model1_ks, mode2_ks], k=2).fit(X, y)
    return obj, X, X_expected

def test_pd(data):
    obj, X, X_expected = data
    X_new = obj.transform(X)
    assert_frame_equal(X_new, X_expected)

@pytest.mark.koalas
def test_ks(data_ks):
    obj, X, X_expected = data_ks
    X_new = obj.transform(X)
    X_new = X_new.to_pandas()
    assert_frame_equal(X_new, X_expected)

def test_pd_np(data):
    obj, X, X_expected = data
    X_numpy_new = obj.transform_numpy(X.to_numpy())
    X_new = pd.DataFrame(X_numpy_new)
    assert_frame_equal(X_new, pd.DataFrame(X_expected.to_numpy()))

@pytest.mark.koalas
def test_ks_np(data_ks):
    obj, X, X_expected = data_ks
    X_numpy_new = obj.transform_numpy(X.to_numpy())
    X_new = pd.DataFrame(X_numpy_new, columns=X_expected.columns)
    assert_frame_equal(X_new, X_expected.astype(np.float64))

def test_init():
    with pytest.raises(TypeError):
        _ = SelectFromModels(models=0, k="a")

    with pytest.raises(TypeError):
        _ = SelectFromModels(models=[XGBClassifier()], k="a")

    class Model:
        pass

    with pytest.raises(TypeError):
        _ = SelectFromModels(models=[Model()], k=2)

max_stars_repo_name: Aditya-Kapadiya/gators

max_issues_repo_name: Aditya-Kapadiya/gators

max_forks_repo_name: Aditya-Kapadiya/gators


text:
import sys

import pytest

@pytest.mark.xfail
def test_simulation_34():
    ...

@pytest.mark.xfail(
    sys.platform.startswith("win"),
    reason="flaky on Windows #42",
    strict=False,
)
def test_login_dialog():
    ...

def check_credentials(name):
    raise NotImplementedError

@pytest.mark.xfail(
    raises=NotImplementedError, reason="will be implemented in #987"
)
def test_credential_check():
    check_credentials("User")  # not implemented yet

class Particle:
    pass

def collide(p1, p2):
    pass

@pytest.mark.xfail(
    run=False, reason="undefined particles cause a crash #625"
)
def test_undefined_particle_collision_crash():
    collide(Particle(), Particle())

def initialize_physics():
    pass

def test_particle_splitting():
    initialize_physics()
    import numpy

    if numpy.__version__ < "1.13":
        pytest.xfail("split computation fails with numpy < 1.13")
    ...


text:
# -*- coding: utf-8 -*-
# <nbformat>3.0</nbformat>

# <codecell>

import numpy as np

# <codecell>

def SetInitialConditions(ICs, ICset = 'Howell', ICtestcase = 0, numPoints = 2000):

    # user inputs
    #ICset = 'Howell'  # 'Sharp' 'Howell' 'Barbee'
    #ICtestcase = 0
    #numPoints = 2000

    # assign simulation variables using specified elements from dictionary of IC's
    mu = ICs[ICset]['mu'][ICtestcase]

    timespan = np.linspace(0, ICs[ICset]['T'][ICtestcase], numPoints)

    initialstate1 = [ICs[ICset]['X'][ICtestcase], 0,                              ICs[ICset]['Z'][ICtestcase],
                     0,                           ICs[ICset]['Ydot'][ICtestcase], 0]

    return mu, timespan, initialstate1

# <codecell>

def InputDataDictionary():

    # create a dictionary for the initial conditions
    ICs = dict()

    # From Sharp, A Collection of Restricted Three-Body Test Problems
    # For problems 1 to 15, mu = 0.012277471 and for problems 16 to 20, mu = 0.000953875
    ICs['Sharp'] = {'mu': np.ones(20),
                    'X':  np.zeros(20),
                    'Z':  np.zeros(20),
                    'Ydot': np.zeros(20),
                    'T':    np.zeros(20)}

    Sharp_X_Z_Ydot_T = np.matrix([[0.994000E+00, 0.0, -0.21138987966945026683E+01, 0.54367954392601899690E+01],
                                  [0.994000E+00, 0.0, -0.20317326295573368357E+01, 0.11124340337266085135E+02],
                                  [0.994000E+00, 0.0, -0.20015851063790825224E+01, 0.17065216560157962559E+02],
                                  [0.997000E+00, 0.0, -0.16251217072210773125E+01, 0.22929723423442969481E+02],
                                  [0.879962E+00, 0.0, -0.66647197988564140807E+00, 0.63006757422352314657E+01],
                                  [0.879962E+00, 0.0, -0.43965281709207999128E+00, 0.12729711861022426544E+02],
                                  [0.879962E+00, 0.0, -0.38089067106386964470E+00, 0.19138746281183026809E+02],
                                  [0.997000E+00, 0.0, -0.18445010489730401177E+01, 0.12353901248612092736E+02],
                                  [0.100000E+01, 0.0, -0.16018768253456252603E+01, 0.12294387796695023304E+02],
                                  [0.100300E+01, 0.0, -0.14465123738451062297E+01, 0.12267904265603897140E+02],
                                  [0.120000E+01, 0.0, -0.71407169828407848921E+00, 0.18337451820715063383E+02],
                                  [0.120000E+01, 0.0, -0.67985320356540547720E+00, 0.30753758552146029263E+02],
                                  [0.120000E+01, 0.0, -0.67153130632829144331E+00, 0.43214375227857454128E+02],
                                  [0.120000E+01, 0.0, -0.66998291305226832207E+00, 0.55672334134347612727E+02],
                                  [0.120000E+01, 0.0, -0.66975741517271092087E+00, 0.68127906604713772763E+02],
                                  [-0.102745E+01, 0.0, 0.40334488290490413053E-01, 0.18371316400018903965E+03],
                                  [-0.976680E+00, 0.0, -0.61191623926410837000E-01, 0.17733241131524483004E+03],
                                  [-0.766650E+00, 0.0, -0.51230158665978820282E+00, 0.17660722897242937108E+03],
                                  [-0.109137E+01, 0.0, 0.14301959822238380020E+00, 0.82949461922342093092E+02],
                                  [-0.110137E+01, 0.0, 0.15354250908611454510E+00, 0.60952121909407746612E+02]])

    ICs['Sharp']['mu'][0:15] *= 0.012277471
    ICs['Sharp']['mu'][15:20] *= 0.000953875
    ICs['Sharp']['X']    = np.array(Sharp_X_Z_Ydot_T[:,0])
    ICs['Sharp']['Z']    = np.array(Sharp_X_Z_Ydot_T[:,1])
    ICs['Sharp']['Ydot'] = np.array(Sharp_X_Z_Ydot_T[:,2])
    ICs['Sharp']['T']    = np.array(Sharp_X_Z_Ydot_T[:,3])

    # From Howell, Three-Dimensional, Periodic, 'Halo' Orbits
    ICs['Howell'] = {'mu':   [0.04, 0.04],
                     'X':    [0.723268,     0.723268],
                     'Z':    [0.040000,    -0.040000],
                     'Ydot': [0.198019,     0.198019],
                     'T':    [1.300177*2.0, 1.300177*2.0]}

    # From Barbee, Notional Mission 4 (Earth-Moon)
    ICs['Barbee'] = {'mu':   [0.012277471],
                     'X':    [0.862307159058101],
                     'Z':    [0.0],
                     'Ydot': [-0.187079489569182],
                     'T':    [2.79101343456226]}  # 12.135 days

    return ICs

# <codecell>

max_stars_repo_name: aerosara/thesis

max_issues_repo_name: aerosara/thesis

max_forks_repo_name: aerosara/thesis


text:
#!/Users/RobbyMiller/opt/anaconda3/bin/python
from __future__ import division

import os
import time
import matplotlib.pyplot as plt
import numpy as np
import copy
import colorsys
import matplotlib.image as mpimg
from IPython.display import display
from IPython.display import Image
import scipy.cluster.hierarchy as sch
from scipy.stats import scoreatpercentile
import scaTools as sca
import mpld3
import pickle
#import cPickle as pickle
from optparse import OptionParser

if not os.path.exists('Outputs/'): os.makedirs('Outputs/')

'''
Summary This script describes the basic flow of the analytical steps in SCA6.0, using the DHFR family as an example.
Here we compare results of the analysis for two different alignments: a PFAM alignment (PFAM PF00186)
and an independent manually curated alignment constructed using a custom database of orthologous sequences (DHFR_PEPM3.an).
Despite differences in the construction, sequence distribution and size of the two alignments,
the sector definition is remarkably consistent: in both cases we arrive at a single sector assembled from six independent components. <br>
For this tutorial, the core calculation scripts should be run from the command line as follows:
>> ./annotate_MSA.py Inputs/PF00186_full.txt -o Outputs/PF00186_full.an -a 'pfam'
>> ./scaProcessMSA.py Inputs/PF00186_full.an -s 1RX2 -c A -f 'Escherichia coli' -t -n
>> ./scaCore.py Outputs/PF00186_full.db
>> ./scaSectorID.py Outputs/PF00186_full.db
>> ./annotate_MSA.py Inputs/DHFR_PEPM3.fasta -o Outputs/DHFR_PEPM3.an -a 'ncbi' -g Inputs/DHFR_PEPM3.gis
>> ./scaProcessMSA.py Inputs/DHFR_PEPM3.an -s 1RX2 -c A -t -n
>> ./scaCore.py Outputs/DHFR_PEPM3.db
>> ./scaSectorID.py Outputs/DHFR_PEPM3.db

Note that we supply annotated alignments for all tutorial scripts (the annotate_pfMSA step is slow, and should only be run once).
O.Rivoire, K.Reynolds and R.Ranganathan 9/2014
'''

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

Dseq = list(); Dsca = list(); Dsect = list()
db = pickle.load(open('Outputs/PF00186_full.db','rb')) # RM! Where you actually import your data.
Dseq.append(db['sequence'])
Dsca.append(db['sca'])
Dsect.append(db['sector'])
db2 = pickle.load(open('Outputs/DHFR_PEPM3.db', 'rb'))
Dseq.append(db2['sequence'])
Dsca.append(db2['sca'])
Dsect.append(db2['sector'])
N_alg = 2
AlgName = ['PFAM', 'Manual']

'''
Read in the results of the above three scripts (scaProcessMSA, scaCore and scaSectorID),
stored as dictionaries in the databases PF00186_full.db and DHFR_PEPM3.db.
To see what variables are stored in each dictionary, use:
	>> print dict.keys()
'''

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

'''

I. Statistical Structure of the Multiple Sequence Alignment (MSA)
We start with a rough comparison of the sequence composition of the two alignments.
Plot a histogram of all pairwise sequence identities (left panel) and a global view of the sequence similarity matrix
(defined by $S\equiv \frac{1}{L}XX^\top$) (right panel).
The PFAM alignment is show in the top row and the manual alignment is shown in the bottom row.
The manual alignment is smaller (644 seqs vs 2000 for PFAM), but both alignments are well-described by a nearly
homogeneous distribution of sequence identities with a mean value of about 35%.
'''

ix = 1
plt.rcParams['figure.figsize'] = 9, 15
for k in range(N_alg):
    # List all elements above the diagonal (i<j):
    listS = [Dsca[k]['simMat'][i,j] for i in range(Dsca[k]['simMat'].shape[0]) \
             for j in range(i+1, Dsca[k]['simMat'].shape[1])]
    #Cluster the sequence similarity matrix
    Z = sch.linkage(Dsca[k]['simMat'],method = 'complete', metric = 'cityblock')
    R = sch.dendrogram(Z,no_plot = True)
    ind = [*map(int, R['ivl'])]
    #Plotting
    plt.rcParams['figure.figsize'] = 14, 4
    plt.subplot(2,2,ix)
    ix += 1
    # listS=np.array(listS)
    # print([[i] for i in listS[0:10]],Dseq[k]['Npos'])
    plt.hist(x=listS, label=Dseq[k]['Npos']/2)
    plt.xlabel('Pairwise sequence identities', fontsize=14)
    plt.ylabel('Number', fontsize=14)
    plt.subplot(2,2,ix)
    ix += 1
    plt.imshow(Dsca[k]['simMat'][np.ix_(ind,ind)], vmin=0, vmax=1)
    plt.colorbar()

plt.show()

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

'''
To examine: (1) sampling differences between the alignments and (2) the relationship between divergence in
sequence similarity and phylogeny for both alignments, we plot the top independent components of the sequence correlation matrix
(after sequence weighting), colored by phylogenetic group. We start by constructing a dictionary of phylogenetic annotations
and checking the representation of sequences in the top taxonomic levels. The annotations are parsed from the sequence headers:
'''

for k in range(N_alg):
    print ("Alignment: " + AlgName[k])
    #construct a dictionary of phylogenetic groups
    annot = dict()
    for i, h in enumerate(Dseq[k]['hd']):
        hs = h.split('|')
        if (len(hs) == 4):
            annot[hs[0]] = sca.Annot(hs[1], hs[2], hs[3].replace('.',''))
        elif (len(hs) == 3):
            annot[hs[0]] = sca.Annot(hs[0],hs[1], hs[2].replace('.',''))

    # Most frequent taxonomic groups:
    atleast = 10
    for level in range(4):
        descr_list = [a.taxo.split(',')[level] for a in annot.values() \
                      if len(a.taxo.split(',')) > level]
        descr_dict = {k:descr_list.count(k) for k in descr_list \
                      if descr_list.count(k)>=atleast}
        print ('\n Level %i:' % level)
        print (descr_dict)

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

'''
Based on this, we select taxonomic groups and colors for representation.
Here, we just start by choosing the broadly well-represented groups.
To see a complete color-coding legend, use:
	>>> sca.figColors()
'''

phylo = list();
fam_names = ['Eukaryota', 'Bacteroidetes', 'Firmicutes', \
             'Actinobacteria', 'Proteobacteria']
col = (0, 0.18, 0.38, 0.6, 0.8)
#Eukaryota = red, Bacteriodetes = yellow, Firmicutes = green,
#Actinobacteria = blue, Proteobacteria = purple
for a in range(N_alg):
    phylo_alg = list()
    for i,k in enumerate(fam_names):
        sf = sca.Unit()
        sf.name = fam_names[i].lower()
        sf.col = col[i]
        sf.items = [j for j,q in enumerate(Dseq[a]['hd'])  \
                    if sf.name in q.lower()]
        phylo_alg.append(sf)
    phylo.append(phylo_alg)

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

'''
Plot the top six independent components of the sequence correlation matrix (with sequence weights);
color-coded by phylogenetic annotation. We compare the phylogenetic sampling for the PFAM alignment (top row)
and manual alignment(bottom row). The data show some very clear seqeunce distinctions based on phylogeny,
and the two alignments seem to differ somewhat in the sequence divergence captured.
In particular, the eukaryotic sequences (in red) seem to form a more distinct group in the manual alignment than in the PFAM alignment.
For the PFAM alignment, the bacteriodetes (yellow) diverge along $U_1$, the actinobacteria (blue) along $U_3$,
the firmicutes (green) along $U_4$ and $U_5$, and a subset of proteobacteria (purple) along $U_6$. For the manual alignment,
the eukaryotes (red) diverge along $U_2$ and $U_6$, the actinobacteria (blue) along $U_4$, the firmicutes (green) along $U_3$,
and a subset of proteobacteria (purple) along $U_5$
'''

plt.rcParams['figure.figsize'] = 9, 8
ix = 1;
for a in range(N_alg):
    U = Dsca[a]['Uica'][1]
    pairs = [[2*i,2*i+1] for i in range(3)]
    for k,[k1,k2] in enumerate(pairs):
        plt.subplot(2,3,ix)
        ix += 1
        sca.figUnits(U[:,k1], U[:,k2], phylo[a])
        #sca.figUnits(U[:,k1], U[:,k2], subfam)
        plt.xlabel(r"${U'}^{(2)}_{%i}$"%(k1+1), fontsize=16)
        plt.ylabel(r"${U'}^{(2)}_{%i}$"%(k2+1), fontsize=16)
    plt.tight_layout()
plt.show()

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

'''
Plot the eigenspectrum of (1) the SCA positional coevolution matrix ($\tilde{C_{ij}}$) (black bars) and
(2) 10 trials of matrix randomization for comparison. This graph is used to choose the number of significant eigenmodes.
Again, we plot the PFAM alignment in the top row and manual alignment in the bottom row for comparison.
Overall the two eigenspectra are remarkably similar: due to small differences in the signficance cutoff,
we define 6 signficant eigenmodes for the PFAM alignment, and 7 for the manual alignment.
'''

for a in range(N_alg):
    plt.rcParams['figure.figsize'] = 9, 4
    hist0, bins = np.histogram(Dsca[a]['Lrand'].flatten(), bins=Dseq[a]['Npos'], \
                               range=(0,Dsect[a]['Lsca'].max()))
    hist1, bins = np.histogram(Dsect[a]['Lsca'], bins=Dseq[a]['Npos'],\
                               range=(0,Dsect[a]['Lsca'].max()))
    plt.subplot(2,1,a+1)
    plt.bar(bins[:-1], hist1, np.diff(bins),color='k')
    plt.plot(bins[:-1], hist0/Dsca[a]['Ntrials'], 'r', linewidth=3)
    plt.tick_params(labelsize=11)
    plt.xlabel('Eigenvalues', fontsize=18); plt.ylabel('Numbers', fontsize=18);
    print ('Number of eigenmodes to keep is %i' %(Dsect[a]['kpos']))
plt.show()
# mpld3.display() #RM! commented out for .py script

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

'''
To define the positions with significant contributions each of the independent components (ICs),
we make a empirical fit for each IC to the t-distribution and select positions with greater than a specified cutoff on the CDF.
We choose $p=0.95$ as our cutoff. Note that since some positions might contribute significantly to more than one IC
(and indication of non-independence of ICs), we apply a simple algorithm to assign such positions to one IC.
Specifically, we assign positions to the IC with which it has the greatest degree of co-evolution.
For brevity, we don't plot the IC fits below (though we do in the other tutorial notebooks),
but do print the list of positions associated with each IC for both alignments.
Comparing between alignments, we can already see some distinctions in the residue positions associated to each component:
IC1 is expanded for the manual alignment, ICs2,4+5 are similar for both alignments, and ICs 3+6 are swapped between the two alignments.
'''

plt.rcParams['figure.figsize'] = 20, 5
for a in range(N_alg):
    print ("alignment: "+AlgName[a])
    for n,ipos in enumerate(Dsect[a]['ics']):
        sort_ipos = sorted(ipos.items)
        ats_ipos = ([Dseq[a]['ats'][s] for s in sort_ipos])
        ic_pymol = ('+'.join(ats_ipos))
        print('IC %i is composed of %i positions:' % (n+1,len(ats_ipos)))
        print(ic_pymol + "\n")

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

'''

To define protein sectors, we examine the structure of the SCA positional correlation matrix with positions contributing
to the top independent components (ICs) ordered by weight (left panel). Again we compare the results between the PFAM alignment (top)
and manual alignment (bottom). This provides a basis to determine/interpret which ICs are truly statistically independent
(defining an independent sector) and which represent hierarchical breakdowns of one sector.
For both alignments, it seems that the ICs reflect a hierarchical break-down of a single sector,
as determined by the high degree of co-evolution in the off-diagonal components (see the dendrogram that follows).
In the right panels the ICs are combined and re-ordered by their contribution to $V_1^p$ to better see this.
'''

sectors = list()
ix = 1
for a in range(N_alg):
    #plot the SCA positional correlation matrix, ordered by contribution
    #to the top ICs
    plt.rcParams['figure.figsize'] = 9, 9
    plt.subplot(2,2,ix); ix +=1;
    plt.imshow(Dsca[a]['Csca'][np.ix_(Dsect[a]['sortedpos'],\
                Dsect[a]['sortedpos'])],vmin=0, vmax=2,\
               interpolation='none',aspect='equal',\
               extent=[0,sum(Dsect[a]['icsize']),0,\
                       sum(Dsect[a]['icsize'])])
    line_index=0
    for i in range(Dsect[a]['kpos']):
        plt.plot([line_index+Dsect[a]['icsize'][i],\
                  line_index+Dsect[a]['icsize'][i]],\
                 [0,sum(Dsect[a]['icsize'])],'w', linewidth = 2)
        plt.plot([0,sum(Dsect[a]['icsize'])],[sum(Dsect[a]['icsize'])\
                        -line_index,sum(Dsect[a]['icsize'])-line_index],\
                 'w', linewidth = 2)
        line_index += Dsect[a]['icsize'][i]

    #combine all the ICs into a single sector and re-sort
    sec_groups = ([k for k in range(Dsect[a]['kpos'])])
    sectors_alg = list()
    s = sca.Unit()
    all_items = list()
    all_Vp = list()
    for i in range(Dsect[a]['kpos']):
        all_items = all_items+Dsect[a]['ics'][i].items
        tmp1 = Dsect[a]['Vpica'][Dsect[a]['ics'][i].items,:]
        all_Vp = all_Vp + list(tmp1[:,0].T)
    svals = list(np.argsort(all_Vp)); svals.reverse()
    s.items = [all_items[i] for i in svals]
    s.col = (1/len(sec_groups))*n
    sectors_alg.append(s)
    sectors.append(sectors_alg)

    #plot the re-ordered matrix
    sortpos = list()
    for s in sectors[a]:
        sortpos.extend(s.items)
    plt.subplot(2,2,ix); ix += 1;
    line_index=0
    plt.imshow(Dsca[a]['Csca'][np.ix_(sortpos, sortpos)], \
               vmin=0, vmax=2,interpolation='none',aspect='equal',\
           extent=[0,len(sortpos),0,len(sortpos)])
    for s in sectors[a]:
        plt.plot([line_index+len(s.items),line_index+len(s.items)],\
                 [0,len(sortpos)],'w', linewidth = 2)
        plt.plot([0,sum(Dsect[a]['icsize'])],[len(sortpos)-line_index, \
                        len(sortpos)-line_index],'w', linewidth = 2)
        line_index += len(s.items)
    plt.tight_layout()
plt.show()

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

'''
The below dendrogram diagrams the relationship between independent components for the PFAM alignment
(the tree for the manual alignment is similar). In this plot, solid lines represent physically contiguous structural units,
and dashed lines indicate spatially fragmented groups of residues.
'''

# i = Image(filename='figs/DHFR_sec_hier.png'); show(i) # RM! For .ipynb files

# RM! Alternative approach for viewing a .png file in .py file.

from PIL import Image
image = Image.open('figs/DHFR_sec_hier.png')
image.show()

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

'''
Print the sector positions (as defined for each alignment), in a format suitable for pyMol, and create
two pyMol sessions with the sectors (and decomposition into independent components) as seperate objects.
'''

for a in range(N_alg):
    print ("Alignment: " + AlgName[a])
    for i,k in enumerate(sectors[a]):
        sort_ipos = sorted(k.items)
        ats_ipos = ([Dseq[a]['ats'][s] for s in sort_ipos])
        ic_pymol = ('+'.join(ats_ipos))
        print('Sector %i is composed of %i positions:' % (i+1,len(ats_ipos)))
        print(ic_pymol + "\n")

sca.writePymol('1RX2', sectors[0], Dsect[0]['ics'], Dseq[0]['ats'],\
               'Outputs/PF00186.pml','A', '../Inputs/', 0)
sca.writePymol('1RX2', sectors[1], Dsect[1]['ics'], Dseq[1]['ats'],\
               'Outputs/DHFR_PEPM3.pml','A', '../Inputs/', 0)

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

'''
As is evident from the position lists above, and as shown below, the structural pattern of the two sectors
and their associated decomposition into independent components is highly similar when compared between the two alignments.
The main difference is that the sector (and independent components) for the manual alignment systematically includes a
few more residue postions.
'''

# i = Image(filename='figs/DHFR_decompv2.png'); i # Use for .ipynb file

image = Image.open('figs/DHFR_decompv2.png') # Alternative approach to view .png in .py file.
image.show()

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

'''
III. The phylogenetic basis of the sector hierarchy
How does the clear phylogenetic heterogeneity in the MSA influence the sector definitions? Since the sector definitions for
the two alignments above are roughly equivalent, we only consider the larger (PFAM) alignment here.
To address this, we take advantage of mathematical methods for mapping between the space of positional and sequence correlations,
as described in Rivoire et al. Using this mapping, we plot the top $k_{pos}$ ICs as 2-D scatter plots with the
corresponding sequence space divergence.
'''

plt.rcParams['figure.figsize'] = 14, 10
pairs= [[0,1],[2,3],[4,5]]
for n,[k1,k2] in enumerate(pairs):
	print('value n:',n)
	print('value k1:',k1)
	print('value k2:',k2)
	plt.subplot(2,3,n+1)
	# print(Dsect[0]['Vpica'][:,3])
	# print(Dsect[0]['Vpica'][:3])
# 	sca.figUnits(Dsect[0]['Vpica'][:,k1], Dsect[0]['Vpica'][:,k2],Dsect[0]['ics'], dotsize = 6) # RM!
# 	plt.xlabel(r'$V^p_{%i}$' % (k1+1), fontsize=16)
# 	plt.ylabel(r'$V^p_{%i}$' % (k2+1), fontsize=16)
# 	plt.subplot(2,3,n+4)
# 	sca.figUnits(Dsect[0]['Upica'][:,k1], Dsect[0]['Upica'][:,k2],\
#                  phylo[0], dotsize = 6)
# 	plt.xlabel(r'$U^p_{%i}$' % (k1+1), fontsize=16)
# 	plt.ylabel(r'$U^p_{%i}$' % (k2+1), fontsize=16)
# plt.tight_layout()
# plt.show()

reason: import numpy,import scipy,from scipy


text:
# Summarize ReFRACtor radiances for baseline forward model run

import numpy
import h5py
from numpy import ndarray

import matplotlib as mpl
mpl.use('agg')
import pylab
from matplotlib import pyplot
from matplotlib import colors
import matplotlib.ticker as mticker

# Sample indicators
clchc = ["#77B0F1","#81BB4E","#E991A1"]
clchcalp = ["#77B0F156","#81BB4E56","#E991A156"]
smpmn = [0, 1016, 2032]
smpmx = [1016, 2032, 3048]

bnms = [r'O$_2$ A-Band',r'Weak CO$_2$',r'Strong CO$_2$']
rdnms = ['o2','weak_co2','strong_co2']

# Control case radiances and wavelengths
rd_ctl = numpy.zeros((3048,))
unc_ctl = numpy.zeros((3048,))
wvln = numpy.zeros((3048,))

ptyp = 'psfc'
tinfo = 'Baseline'

ctlfl = 'land_fmout_control.h5'
f = h5py.File(ctlfl)
for j in range(3):
    sst0 = smpmn[j]
    sfn0 = smpmx[j]
    vrnm = '/SoundingMeasurements/noiseless_radiance_%s' % (rdnms[j])
    rd_ctl[sst0:sfn0] = f[vrnm][:] * 1.0e-20
    vrnm = '/SoundingMeasurements/radiance_uncert_%s' % (rdnms[j])
    unc_ctl[sst0:sfn0] = f[vrnm][:] * 1.0e-20
    wvnm = '/SoundingMeasurements/wavelength_%s' % (rdnms[j])
    wvln[sst0:sfn0] = f[wvnm][:]
f.close()

fig = pyplot.figure(figsize=(15,5))
for bnd in range(3):
    smpsq = numpy.arange(smpmn[bnd],smpmx[bnd])
    smpidx = smpsq - 1016 * bnd

    p1spt = bnd + 1
    p1 = pyplot.subplot(1,3,p1spt)
    pl1, = p1.plot(wvln[smpsq],rd_ctl[smpsq],ls='-',color=clchc[bnd])
    p1.xaxis.grid(color='#898989',linestyle='dotted')
    p1.yaxis.grid(color='#898989',linestyle='dotted')
    p1.set_xlabel(r'Wavelength [$\mu$m]',size=11)
    p1.set_ylabel('Radiance x 1e20',size=11)
    if bnd == 0:
        p1.set_xticks([0.758,0.760,0.762,0.764,0.766,0.768,0.770])
        p1.set_xticklabels(['','0.760','','0.764','','0.768',''])
    for lb in p1.xaxis.get_ticklabels():
        lb.set_fontsize(11)
    for lb in p1.yaxis.get_ticklabels():
        lb.set_fontsize(11)
    tstr = r'%s OCO-2 Radiance' % (bnms[bnd])
    pyplot.title(tstr)

fig.subplots_adjust(bottom=0.1,top=0.9,left=0.05,right=0.95, \
                    hspace=0.35,wspace=0.35)
fig.suptitle(tinfo,fontsize=14)
fnm = 'Refractor_OCO2_FMBaseline.png'
pyplot.savefig(fnm)
pyplot.close()

reason: import numpy,from numpy


text:
# ************************ General Information ************************
'''
VERSION:
-------

Version (by release date): 2018-07-26

DEVELOPER INFORMATION:
---------------------

Name: Omid Mirzaei
Laboratory: Computer Security Lab (COSEC)
University: Universidad Carlos III de Madrid
Website: https://cosec.inf.uc3m.es/~omid-mirzaei/

PUBLICATION:
-----------

AndrODet: An Adaptive Android Obfuscation Detector
O. Mirzaei, J. M. de Fuentes, J. E. Tapiador, L. Gonzalez-Manzano
Future Generation Computer Systems, Elsevier (January 2019)

COPYRIGHT NOTICE:
----------------

All rights reserved for the above developer and research center.
Please, take a look at the "License.txt" file for more detailed information regarding the usage and distribution of these source codes.

ACKNOWLEDGEMENT:
---------------

This work has been partially supported by the:
MINECO grant TIN2016-79095-C2-2-R (SMOG-DEV);
CAM grant S2013/ICE-3095 (CIBERDINE);
co-funded with European FEDER funds;
partially supported by the UC3M's grant Programa de Ayudas para la Movilidad.
The authors would like to thank the Allatori technical team for its valuable assistance, and, also, the authors of the AMD and PraGuard datasets which made their repositories available to us.
'''
# ************************ End of General Information ************************

# ************************ Module Information ************************
'''
MAIN FUNCTIONALITY:
------------------

This module extracts features from each Android application in the stream, test the system, and, also, trains it incrementally.

ARGUMENTS:
---------

-a:     Directory of Android applications (.apk files)
-d:     Directory of dexdump disassembler.
-g:     Directory of androguard tool.
-o:     Directory of output.

USAGE:
-----

python AndrODet_MOA.py -a '/Directory/of/apps' -d '/Directory/of/dexdump' -g '/Directory/of/androguard' -o '/Directory/of/output'
'''
# ************************ End of Module Information  ************************

# ************************ Importing Modules ************************

import os
from sets import Set
import subprocess
import multiprocessing
from optparse import OptionParser
from feature_extraction import EFI
from feature_extraction import EFS
from feature_extraction import EFC
from sklearn.metrics import confusion_matrix, accuracy_score
import numpy as np

# ************************ End of Importing Modules ************************

# ************************ Initialization ************************

JAVA_PATH = '/Library/Internet Plug-Ins/JavaAppletPlugin.plugin/Contents/Home/bin/java'
Home_Dir = os.path.curdir                               # Home directory
MOA_CP = os.path.join(Home_Dir, 'MOA')                  # MOA directory
n_procs = 20                                            # Number of processes
num_features_IR = 15                                    # Number of features for IR
num_features_SE = 8                                     # Number of features for SE
num_features_CF = 7                                     # Number of features for CF
Real_Classes = {}                                       # Real classes of apps

# --------------- Setting command-line options ---------------
option_1 = { 'name' : ('-a', '--apps_dir'), 'help' : 'Directory of apk files', 'nargs' : 1 }
option_2 = { 'name' : ('-d', '--dexdump_dir'), 'help' : 'Directory of dexdump', 'nargs' : 1 }
option_3 = { 'name' : ('-g', '--androguard_dir'), 'help' : 'Directory of androguard', 'nargs' : 1 }
option_4 = { 'name' : ('-o', '--output_dir'), 'help' : 'Directory of output', 'nargs' : 1 }

options = [option_1, option_2, option_3, option_4]
# --------------- End of Setting command-line options ---------------

# --------------- Configuring MOA settings ---------------

class IR_Detector():

    def __init__(self):
        self.learner_IR = 'meta.LeveragingBag'                          # Learning algorithm for identifier renaming
        self.class_col_num_IR = str(num_features_IR + 1)                # Class column number in .arff file

    def extract_features(self, appfile):
        dirname, filename = os.path.split(appfile)
        dirname = os.path.join(dirname, 'apps_features')
        app_name = filename[:-4]
        features_IDs = EFI.extract_features(appfile, options.apps_dir, options.dexdump_dir, dirname)
        return features_IDs

    def test(self, arff_file):
        subprocess.call([JAVA_PATH, '-cp', os.path.join(MOA_CP, 'moa.jar'), \
                        '-javaagent:' + os.path.join(MOA_CP, 'sizeofag-1.0.0.jar'), 'moa.DoTask', \
                        'EvaluatePrequential', '-l', '(' + self.learner_IR, '-s', '20)', \
                        '-s', '(ArffFileStream', '-f', arff_file, '-c', self.class_col_num_IR + ')', \
                        '-i', '-1', '-f', '1', '-o', os.path.join(options.output_dir, 'predictions_IR_module')])
        output_file = open(os.path.join(options.output_dir, 'predictions_IR_module'), 'rb')
        result = output_file.readlines()
        return result

class SE_Detector():

    def __init__(self):
        self.learner_SE = 'meta.LeveragingBag'                           # Learning algorithm for string encryption
        self.class_col_num_SE = str(num_features_SE + 1)                 # Class column number in .arff file

    def extract_features(self, appfile):
        dirname, filename = os.path.split(appfile)
        dirname = os.path.join(dirname, 'apps_features')
        app_name = filename[:-4]
        features_STs = EFS.extract_features(appfile, options.apps_dir, options.dexdump_dir, dirname)
        return features_STs

    def test(self, arff_file):
        subprocess.call([JAVA_PATH, '-cp', os.path.join(MOA_CP, 'moa.jar'), \
                        '-javaagent:' + os.path.join(MOA_CP, 'sizeofag-1.0.0.jar'), 'moa.DoTask', \
                        'EvaluatePrequential', '-l', '(' + self.learner_SE, '-s', '20)', \
                        '-s', '(ArffFileStream', '-f', arff_file, '-c', self.class_col_num_SE + ')', \
                        '-i', '-1', '-f', '1', '-o', os.path.join(options.output_dir, 'predictions_SE_module')])
        output_file = open(os.path.join(options.output_dir, 'predictions_SE_module'), 'rb')
        result = output_file.readlines()
        return result

class CF_Detector():

    def __init__(self):
        self.learner_CF = 'meta.LeveragingBag'                            # Learning algorithm for control flow obfuscation
        self.class_col_num_CF = str(num_features_CF + 1)                  # Class column number in .arff file

    def extract_features(self, appfile):
        dirname, filename = os.path.split(appfile)
        dirname = os.path.join(dirname, 'apps_features')
        app_name = filename[:-4]
        features_CFs = EFC.extract_features(appfile, options.apps_dir, options.androguard_dir, options.dexdump_dir, dirname)
        return features_CFs

    def test(self, arff_file):
        subprocess.call([JAVA_PATH, '-cp', os.path.join(MOA_CP, 'moa.jar'), \
                        '-javaagent:' + os.path.join(MOA_CP, 'sizeofag-1.0.0.jar'), 'moa.DoTask', \
                        'EvaluatePrequential', '-l', '(' + self.learner_CF, '-s', '20)', \
                        '-s', '(ArffFileStream', '-f', arff_file, '-c', self.class_col_num_CF + ')', \
                        '-i', '-1', '-f', '1', '-o', os.path.join(options.output_dir, 'predictions_CF_module')])
        output_file = open(os.path.join(options.output_dir, 'predictions_CF_module'), 'rb')
        result = output_file.readlines()
        return result

# --------------- End of Configuring MOA settings ---------------

# ********************* End of Initialization *********************

# ********************* Functions *********************

def set_learners():
    # --------------- Setting learner parameters ---------------

    IR_module = IR_Detector()
    SE_module = SE_Detector()
    CF_module = CF_Detector()

    # --------------- End of Setting learner parameters ---------------

    return IR_module, SE_module, CF_module

def feature_extraction(appfile, IR_module, SE_module, CF_module):
    dirname, filename = os.path.split(appfile)
    dirname = os.path.join(dirname, 'apps_features')
    app_name = filename[:-4]
    try:
        # --------------- Extracting features ---------------

        features_IR = IR_module.extract_features(appfile)
        features_SE = SE_module.extract_features(appfile)
        features_CF = CF_module.extract_features(appfile)

        # --------------- End of Extracting features ---------------

        if features_IR and len(features_IR) == num_features_IR and features_SE and len(features_SE) == num_features_SE and features_CF and len(features_CF) == num_features_CF:
            features_IR.append(Real_Classes[app_name][0])
            features_SE.append(Real_Classes[app_name][1])
            features_CF.append(Real_Classes[app_name][2])
            return app_name, features_IR, features_SE, features_CF
    except:
        print 'features extraction failed for app', appfile

def detect_obfuscation(IR_arff_fie, SE_arff_fie, CF_arff_fie):
    try:
        # --------------- Testing the learner ---------------

        predict_output_IR = IR_module.test(IR_arff_fie)
        predict_output_SE = SE_module.test(SE_arff_fie)
        predict_output_CF = CF_module.test(CF_arff_fie)

        # --------------- End of Testing the learner ---------------

        return predict_output_IR, predict_output_SE, predict_output_CF
    except:
        print 'System could not successfully analyze app %s!' %app_name
        pass

def confusion_matrix_update(conf_matrix, real_classes, predicted_classes):

    row_idx = int(real_classes, 2)
    col_idx = int(predicted_classes, 2)

    conf_matrix[row_idx][col_idx] += 1
    return conf_matrix

# ********************* End of Functions *********************

# ********************* Main Body *********************

if __name__ == '__main__':

    parser = OptionParser()
    for option in options:
        param = option['name']
        del option['name']
        parser.add_option(*param, **option)
    options, arguments = parser.parse_args()

    if not os.path.exists(options.output_dir):
        os.mkdir(options.output_dir)

    IR_module, SE_module, CF_module = set_learners()

    all_apks = Set([])                  # To discard possible redundant apk files
    for root, directories, filenames in os.walk(options.apps_dir):
        for filename in filenames:
            if '.DS_Store' not in filename  and '.apk' in filename:
                all_apks.add(os.path.join(root,filename))
                if filename[:-4] not in Real_Classes.keys():
                    Real_Classes[filename[:-4]] = np.zeros(3, dtype=int)
                if '/IR/YES' in root:
                    Real_Classes[filename[:-4]][0] = 1
                if '/SE/YES' in root:
                    Real_Classes[filename[:-4]][1] = 1
                if '/CF/YES' in root:
                    Real_Classes[filename[:-4]][2] = 1

    pool = multiprocessing.Pool(n_procs)
    results = [pool.apply_async(feature_extraction, [appfile, IR_module, SE_module, CF_module]) for appfile in all_apks]
    pool.close()
    pool.join()

    features_IR = []
    features_SE = []
    features_CF = []
    processed_apps = []
    for res in results:
        try:
            if len(res.get()) == 4:
                app_name = res.get()[0]
                features_IR.append(res.get()[1])
                features_SE.append(res.get()[2])
                features_CF.append(res.get()[3])

                processed_apps.append(app_name)
        except:
            pass

    features_file_IR_arff = open(os.path.join(options.output_dir, 'features_IR.arff'), 'wb')
    features_file_SE_arff = open(os.path.join(options.output_dir, 'features_SE.arff'), 'wb')
    features_file_CF_arff = open(os.path.join(options.output_dir, 'features_CF.arff'), 'wb')

    EFI.save_features_to_arff(features_IR, features_file_IR_arff)
    EFS.save_features_to_arff(features_SE, features_file_SE_arff)
    EFC.save_features_to_arff(features_CF, features_file_CF_arff)

    features_file_IR_arff.close()
    features_file_SE_arff.close()
    features_file_CF_arff.close()

    predict_output_IR, predict_output_SE, predict_output_CF =  detect_obfuscation(os.path.join(options.output_dir, 'features_IR.arff'), \
                                                                                os.path.join(options.output_dir, 'features_SE.arff'), \
                                                                                os.path.join(options.output_dir, 'features_CF.arff'))

    conf_matrix = np.zeros((8, 8), dtype=int)
    for idx in range(0, len(predict_output_IR)):
        real_classes = str(Real_Classes[processed_apps[idx]][0]) + str(Real_Classes[processed_apps[idx]][1]) + str(Real_Classes[processed_apps[idx]][2])
        predicted_classes = predict_output_IR[idx].split(',')[0] + predict_output_SE[idx].split(',')[0] + predict_output_CF[idx].split(',')[0]
        conf_matrix = confusion_matrix_update(conf_matrix, real_classes, predicted_classes)

    print 'Confusion Matrix:'
    print conf_matrix
# ********************* End of Main Body *********************

max_stars_repo_name: rpanah/AndrODet

max_issues_repo_name: rpanah/AndrODet

max_forks_repo_name: rpanah/AndrODet


text:
import numpy as np
import tensorflow as tf

def loadData():
    # Loading my data
    trainData = np.load('data2D.npy')

    # random shuffle trainData and trainTarget
    np.random.seed(521)
    np.random.shuffle(trainData)
    validData_random = trainData[:3333, :]
    trainData_random = trainData[3334:, :]
    return trainData_random, validData_random

def LOSS(X,Mu):
    # compute the loss function
    #X:trainData [10000 2]
    #Mu: center value matrix [3 2]
    dist = tf.reduce_sum(tf.square(X[:, tf.newaxis] - Mu), 2)
    loss = tf.reduce_sum(tf.reduce_min(dist,1))
    return loss, dist

def evaluate(dist,K):
    percentage = []
    _, indices = tf.nn.top_k(-dist)
    indices = tf.reshape(indices,[-1])
    for i in xrange(K):
        data = tf.cast(i*tf.ones((6666,)),tf.int32)
        correct = tf.equal(indices,data)
        percentage.append(tf.reduce_sum(tf.cast(correct, tf.float64)) / 6666)

    return percentage

def buildGraph(learingRate,K):
    # Variable creation
    Mu = tf.Variable(tf.random_normal(shape=[K,2],stddev=1.0,dtype=tf.float64),name = 'Mu')
    X = tf.placeholder(dtype=tf.float64, shape=[6666, 2], name='trainData')

    # Graph and Error definition
    totalLoss, dist = LOSS(X,Mu)

    # Training mechanism
    optimizer = tf.train.AdamOptimizer(learning_rate = learingRate,beta1=0.9, beta2=0.99,epsilon=1e-5)
    train = optimizer.minimize(loss=totalLoss)
    return X, Mu, totalLoss, train, dist

def runKmeans(learningRate, K, epoch):
    # Build computation graph
    X, Mu, totalLoss, train, dist = buildGraph(learningRate,K)

    # Loading my data
    trainData, validData = loadData()

    # Initialize session
    init = tf.global_variables_initializer()

    sess = tf.InteractiveSession()
    sess.run(init)

    initialMu = sess.run(Mu)

    # set for logging
    j = 0

    # Training model
    logging = np.zeros((epoch, 2))

    for step in xrange(0, epoch):
        # trainprocess on training set
        _,currentDist, currentMu, trainLoss = sess.run([train, dist, Mu, totalLoss],feed_dict={X: trainData})

        validLoss, _ =sess.run(LOSS(validData,currentMu))
        logging[j] = [trainLoss,validLoss]
        j += 1

    print ('Final valid loss:',validLoss,'when K=',K)

    return logging

max_stars_repo_name: ouyangyike/Inference-Algorithm

max_issues_repo_name: ouyangyike/Inference-Algorithm

max_forks_repo_name: ouyangyike/Inference-Algorithm


text:
import sys, re
import pandas as pd
import numpy as np

import torch
from torch.utils.data import DataLoader
from torchvision import models, transforms
from sklearn.metrics import roc_auc_score

basepath = '/USER/INFERENCE/CANCER/'
sys.path.append(basepath + 'src/')
from dataloader_cancer import ImageDataset
from utils import seed_everything, my_model

seed = 20
seed_everything(seed)

def infer(model, data_loader, device):
    codes = []
    preds = np.array([])

    print_once = 0
    model.eval()
    with torch.no_grad():
        for i,data in enumerate(data_loader, 0):
            inputs, code = data['image'], data['code']
            if not print_once:
                print('** image size: ', inputs.shape)
                print_once = 1
            inputs = inputs.to(device)

            outputs = model.forward(inputs)
            proba = torch.nn.Softmax(dim=1)(outputs)

            if i%100 == 99:
                print(i+1)

            codes.extend(code)
            preds = np.hstack((preds, np.array(proba.cpu())[:,1]))

    return codes, preds

def get_transform(option):
    normalize = transforms.Normalize(mean=[0.485, 0.456, 0.406], std=[0.229, 0.224, 0.225])
    if option == 0:
        transformer = transforms.Compose([transforms.ToTensor(), normalize])
    elif option == 1:
        transformer = transforms.Compose([transforms.RandomHorizontalFlip(1), transforms.ToTensor(), normalize])
    elif option == 2:
        transformer = transforms.Compose([transforms.RandomVerticalFlip(1), transforms.ToTensor(), normalize])
    return transformer

def main(filename, savename, model_type):
    print('** Input filenames: ', filename)
    num_model = len(filename)

    if model_type!=None and num_model != len(model_type):
        print('## ERROR: number of filenames and model types are not matched!')
    if model_type == None:
        print('## No model type is defined.')
        print('## --> Set as default (efficientnet-b0)')
        model_type = [0]*num_model

    model_names = ['efficientnet-b0', 'efficientnet-b5']

    print('** Model types: ', model_type)

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

    normalize = transforms.Normalize(mean=[0.485, 0.456, 0.406], std=[0.229, 0.224, 0.225])
    transformer = transforms.Compose([transforms.ToTensor(), normalize])

    averaged = []
    for i in range(len(filename)):
        fn = filename[i]
        mtype = model_names[model_type[i]]
        print(f'\n{i+1}/{num_model}... infer ... {fn}.pth ... {mtype}')
        if mtype == 'efficientnet-b0': small_size = 0
        else: small_size = 1

        test_dataset = ImageDataset(mode='test', transform=transformer, smaller=small_size)
        test_loader = torch.utils.data.DataLoader(test_dataset, batch_size=1, shuffle=False)

        model = my_model(3, 2, mtype)
        model.load_state_dict(torch.load(f'{basepath}weights/{fn}.pth'))
        model.to(device)

        code, pred = infer(model, test_loader, device)

        averaged.append(pred)

        pd.DataFrame({'filename': code, 'pred': pred}).to_csv(f'{basepath}results/{fn}.csv', index=False)

    if len(filename) > 1:
        averaged = np.mean(np.array(averaged), axis=0)

        print(averaged)
        pd.DataFrame({'filename': code, 'pred': averaged}).to_csv(f'{basepath}results/{savename}.csv', index=False)

if __name__ == "__main__":
    import argparse
    parser = argparse.ArgumentParser(description = "Step2_ Write result file")
    parser.add_argument(dest='filenames', type=str, nargs='+', help="Filename for model.pth and result.csv")
    parser.add_argument('-m', dest='model_types', type=int, nargs='+', required=False, help="The type of the model")
    parser.add_argument('-s', dest='savename', type=str, required=False, help="Filename to save in case of ensemble")
    #parser.add_argument('-s', dest='small_size', type=int, nargs='+', required=False, help="Downscaling input")
    args = parser.parse_args()

    main(args.filenames, args.savename, args.model_types)#, args.small_size)


text:
#!/usr/bin/env python
#-----------------------------------------------------------------------
#
# Description:
#
# GDAL based GeoTIFF utilities
#
#-------------------------------------------------------------------------------
#
# Project: EOxServer <http://eoxserver.org>
# Authors: Martin Paces <martin.paces@iguassu.cz>
#
#-------------------------------------------------------------------------------
# Copyright (C) 2011 Iguassu Software Systems, a.s
#
# Permission is hereby granted, free of charge, to any person obtaining a copy
# of this software and associated documentation files (the "Software"), to deal
# in the Software without restriction, including without limitation the rights
# to use, copy, modify, merge, publish, distribute, sublicense, and/or sell
# copies of the Software, and to permit persons to whom the Software is
# furnished to do so, subject to the following conditions:
#
# The above copyright notice and this permission notice shall be included in all
# copies of this Software or works derived from this Software.
#
# THE SOFTWARE IS PROVIDED "AS IS", WITHOUT WARRANTY OF ANY KIND, EXPRESS OR
# IMPLIED, INCLUDING BUT NOT LIMITED TO THE WARRANTIES OF MERCHANTABILITY,
# FITNESS FOR A PARTICULAR PURPOSE AND NONINFRINGEMENT. IN NO EVENT SHALL THE
# AUTHORS OR COPYRIGHT HOLDERS BE LIABLE FOR ANY CLAIM, DAMAGES OR OTHER
# LIABILITY, WHETHER IN AN ACTION OF CONTRACT, TORT OR OTHERWISE, ARISING FROM,
# OUT OF OR IN CONNECTION WITH THE SOFTWARE OR THE USE OR OTHER DEALINGS IN
# THE SOFTWARE.
#-------------------------------------------------------------------------------

import sys
import os.path
from numpy import linalg
from numpy import array
from numpy import dot
from numpy import argmin, argmax

from qhull import chull2D_qhull

from eoxserver.contrib import osr, gdal

#-------------------------------------------------------------------------------

class GeoTransform(object) :

    def __init__( self , gt ) :

        self.xy = array( (gt[0],gt[3]) , 'float64' )
        self.af = array( ((gt[2],gt[1]),(gt[5],gt[4])) , 'float64' )
        self.ab = linalg.inv( self.af )

    def rc2xy( self , rc ) :
        rc = array( rc , 'float64' )
        xy = self.xy + dot( self.af , rc )
        return xy

    def xy2rc( self , xy ) :
        xy = array( xy , 'float64' )
        rc = dot( self.ab , xy - self.xy )
        return rc

#-------------------------------------------------------------------------------

class OSRTransform( object ) :

    def __init__( self , src , dst ) :

        self.fwd = osr.CoordinateTransformation( src , dst )
        self.bwd = osr.CoordinateTransformation( dst , src )

    def src2dst( self , xy ) :
        return array( self.fwd.TransformPoint( xy[0] , xy[1] ) , 'float64' )

    def dst2src( self , xy ) :
        return array( self.bwd.TransformPoint( xy[0] , xy[1] ) , 'float64' )

#-------------------------------------------------------------------------------

class GDalInfo( object ) :
    """ this class holds the GDAL information of the rester file """

    def __init__( self , fname ) :

        # get the dataset and its params
        gdal.AllRegister()
        ds = gdal.Open(fname)

        self.fileName = fname
        self.driverName =  ds.GetDriver().LongName
        self.size       = (ds.RasterYSize,ds.RasterXSize,ds.RasterCount)
        self.GCPCount   = ds.GetGCPCount()
        self.GCPProjection = ds.GetGCPProjection()
        self.Projection = ds.GetProjection()
        self.ProjectionRef = ds.GetProjectionRef()
        self.GeoTransform = ds.GetGeoTransform()
        self.GCP         = ds.GetGCPs()
        self.isRectified  = bool( self.Projection )
        self.isReferenceable = bool( self.GCPProjection) and ( self.GCPCount > 0 )

    def __str__( self ) :

        a = []
        a.append( "GDalInfo:" )
        a.append( "\tFile:         \t %s" % str(self.fileName) )
        a.append( "\tDriver:       \t %s" % str(self.driverName) )
        a.append( "\tSize:         \t %s" % str( self.size ) )
        a.append( "\tNo.GCP        \t %s" % str( self.GCPCount ) )
        a.append( "\tProjectionRef:\t %s" % str( self.ProjectionRef  ) )
        a.append( "\tProjection:   \t %s" % str( self.Projection ) )
        a.append( "\tGeoTransform: \t %s" % str( self.GeoTransform ) )
        a.append( "\tGCP Projection:\t %s" % str( self.GCPProjection ) )
        a.append( "" )

        return "\n".join( a )

#-------------------------------------------------------------------------------

def getFootprint( info , **kvarg ) :
    """ extract geotiff image footprint

        info            - [instance of GDalInfo class] image descriptor
        numberOfPoints  - [integer] number of points per single image size
                          (20 by default, applicable to Rectified images only)
        delimiter       - [string] separator delimiting the footprint coordinates (<space> by default)
        repeatFirst     - [boolean] if True the first point is repeated as last one to close the loop
    """

    if info.isRectified :
        return getFootprintRect( info , **kvarg )
    elif info.isReferenceable :
        return getFootprintRef( info , **kvarg )

    return None

def getFootprintRect( info , numberOfPoints = 20, delimiter = " " , repeatFirst = False ) :
    """ extract geotiff image footprint of a Rectified image

        info            - [instance of GDalInfo class] image descriptor
        numberOfPoints  - [integer] number of points per single image size
                          (20 by default)
        delimiter       - [string] separator delimiting the footprint coordinates (<space> by default)
        repeatFirst     - [boolean] if True the first point is repeated as last one to close the loop
    """

    # -----------------------------
    # projections conversion

    size = ( info.size[0] , info.size[1] )

    gt = GeoTransform( info.GeoTransform )

    cr_src = osr.SpatialReference()
    cr_dst = osr.SpatialReference()

    cr_src.ImportFromWkt( info.Projection )
    cr_dst.SetWellKnownGeogCS( "WGS84" )

    trn = OSRTransform(cr_src,cr_dst)

    # -----------------------------
    # extract footprint (clockwise)

    foot = []

    # top (L2R)
    for i  in xrange( 0 , numberOfPoints ) : foot.append( tuple( trn.src2dst(gt.rc2xy((0,size[1]*float(i)/float(numberOfPoints))))[:2] ) )
    # right (T2B)
    for i  in xrange( 0 , numberOfPoints ) : foot.append( tuple( trn.src2dst(gt.rc2xy((size[0]*float(i)/float(numberOfPoints),size[1])))[:2] ) )
    # bottom (R2L)
    for i  in xrange( 0 , numberOfPoints ) : foot.append( tuple( trn.src2dst(gt.rc2xy((size[0],size[1]*float(numberOfPoints-i)/float(numberOfPoints))))[:2] ) )
    # left (B2T)
    for i  in xrange( 0 , numberOfPoints ) : foot.append( tuple( trn.src2dst(gt.rc2xy((size[0]*float(numberOfPoints-i)/float(numberOfPoints),0)))[:2] ) )

    if ( repeatFirst ) : foot.append( foot[0] )

    # things to be done
    # 1) POLAR PROJECTIONS
    # 2) NON-POLAR FOOTPRINT CROSSING +/-180 MERIDIAN

    return delimiter.join(map( lambda x : "%.6f%s%.6f"%(x[0],delimiter,x[1]) , foot ))

#-------------------------------------------------------------------------------

def getFootprintRef( info , numberOfPoints = 20, delimiter = " " , repeatFirst = False ) :
    """ extract geotiff image footprint of a Referenceableimage

        info            - [instance of GDalInfo class] image descriptor
        numberOfPoints  - [integer] dummy placeholder - not used (20 by default)
        delimiter       - [string] separator delimiting the footprint coordinates (<space> by default)
        repeatFirst     - [boolean] if True the first point is repeated as last one to close the loop
    """
    # NOTE: numberOfPoints parameter is ignored

    D = array( map( lambda gcp : ( gcp.GCPLine , gcp.GCPPixel ) , info.GCP ) )

    H = chull2D_qhull( D , eps = 1e-3 , clock_wise = True )

    if ( repeatFirst ) : H.append( H[0] )

    # things to be done
    # 1) POLAR PROJECTIONS
    # 2) NON-POLAR FOOTPRINT CROSSING +/-180 MERIDIAN

    return delimiter.join(map( lambda i : "%.6f%s%.6f"%(info.GCP[i].GCPX,delimiter,info.GCP[i].GCPY) , H ))

#-------------------------------------------------------------------------------

# simple test
if __name__ == "__main__" :

    fname = os.path.abspath( sys.argv[1] )

    info = GDalInfo( fname )

    #print info
    #print "WGS84 Footprint: "
    print getFootprint( info )

max_stars_repo_name: ESA-VirES/eoxserver

max_issues_repo_name: ESA-VirES/eoxserver

max_forks_repo_name: ESA-VirES/eoxserver

reason: from numpy


text:
#!/usr/bin/env python3
# -*- coding: utf-8 -*-

#===============================================================================
# PROGRAM METADATA
#===============================================================================
__author__ = 'xinwei <xwei97@zju.edu.cn>'
__copyright__ = ''
__license__ = ''
__date__ = 'Tue Mar 23 14:30:10 2021'
__version__ = '0.1'

USAGE = \
"""
    %prog [options] -i [path of input.kreport] -o [path of output.kreport]

    Description:
        Order by the rank of classification, put the rows with same rank together.

    Dependency:
        pandas, os

    Examples:

    AUTHOR: Xin Wei <xwei97@zju.edu.cn>
    SINCE: Tue Mar 23 14:30:10 2021
    VERSION: 0.0.0

    Kraken2 output format:
    1. Percentage of fragments covered by the clade rooted at this taxon
    2. Number of fragments covered by the clade rooted at this taxon
    3. Number of fragments assigned directly to this taxon
    4. A rank code, indicating (U)nclassified, (R)oot, (D)omain, (K)ingdom, (P)hylum, (C)lass, (O)rder, (F)amily, (G)enus, or (S)pecies.
    Taxa that are not at any of these 10 ranks have a rank code that is formed by using the rank code of the closest ancestor rank with a number indicating the distance from that rank.
    E.g., "G2" is a rank code indicating a taxon is between genus and species and the grandparent taxon is at the genus rank.
    5. NCBI taxonomic ID number
    6. Indented scientific name
"""

#===============================================================================
# IMPORT STATEMENTS
#===============================================================================
import pandas as pd
import os
import numpy as np
import re
import sys
import time
from optparse import OptionParser, OptionGroup

#===============================================================================
# MAIN
#===============================================================================
def main(argv = sys.argv[1:]):
    parser = OptionParser( usage=USAGE, version = "Version: " + __version__)
	# Compulsory arguments
    compOptions = OptionGroup(parser, "Compulsory parameters", "There options are compulsory, and may be supplied in any order.")
    compOptions.add_option("-i", "--indir", type = "string", metavar = "DIR", help = "the directory of .kreport file generated by kraken2")
    compOptions.add_option("-o", "--outdir", type = "string", metavar = "DIR", help = "the directory of reordered .kreport file")
    parser.add_option_group(compOptions)
    try:
        options, x = parser.parse_args(argv)
    except:
        parser.error('[FATAL]: parsing failure.\n')
        exit(1)

    if not options.indir:
        parser.error('[FATAL]: input directory is not accessible.\n')
        exit(1)
    if not options.outdir:
        parser.error('[FATAL]: output directory is not accessible.\n')
        exit(1)
    inDir = options.indir
    outDir = options.outdir
    fileList = os.listdir(inDir)
    for i in fileList:
        if ('.report' in i) & ('_sorted.report' not in i):
             df = pd.read_csv(inDir + '/' + i, sep = '\t', header = None, names = ['percentage', 'reads_covered', 'reads_assigned', 'rank', 'taxid', 'name'])
             df['name'] = [s.strip() for s in df['name']]
             uniRank = set(df['rank'])
             uniRank_list = [re.split(r'(\d+)', u, maxsplit=1) for u in uniRank]
             uni = list()
             for j in uniRank_list:
                 if len(j) == 1:
                     t = (j[0], 0)
                 if len(j) == 3:
                     t = (j[0], int(j[1]))
                 uni.append(t)
             d = {'U':1, 'R':2, 'S':3, 'G':4, 'F':5, 'O':6, 'C':7, 'P':8, 'K':9, 'D':10}
             uni.sort(key=lambda x:(d[x[0]], x[1]))
             custom_ord = []
             ten = []
             for x,y in uni:
                 if y == 0:
                     z = x
                     ten.append(z)
                 else:
                     z = x + str(y)
                     custom_ord.append(z)
             custom_ord = ten + custom_ord

             df['rank'] = pd.Categorical(df['rank'], custom_ord)
             df.sort_values(by = ['rank', 'percentage', 'reads_assigned', 'reads_covered', 'taxid'], ascending=[True] + [False]*4, inplace = True)
             df.to_csv(outDir + '/' + re.sub(r'\.report$', '', i) + '_sorted.report', sep='\t', index=False, header=True)
        else:
            pass

if __name__ == '__main__':
    t0 = time.time()
    main()
    sys.stderr.write('Elapsed time to run reorder_kraken2_report.py: {} s\n'.format( (time.time()-t0) ) )


text:
'''Code relating to the mask learning agents use to understand the game'''
import numpy
from abc import ABC, abstractmethod

class Mask(ABC):
    '''Abstract class identifing the functions a mask class needs'''

    def __init__(self, name, boardSize=4, maxTile=15):
        '''Init the mask class
        input:
            name: Name of mask used for tag
            boardSize: Size of the board
            maxTile: log2 of max tile that can appear on board'''
        self.name = name
        self.boardSize = boardSize
        self.maxTile = maxTile

    @abstractmethod
    def getNumTuples(self):
        '''Return the number of tuples used to describe each state'''
        pass

    @abstractmethod
    def getMaxTupleNum(self):
        '''Return the largest number corresponding to a tuple'''
        pass

    @abstractmethod
    def getTupleNums(self, state):
        '''Transforms a state into its tuple number representation
        input:
            state: state to transform
        output: array of tuple nums corresponding to state'''
        pass

    def getTag(self):
        '''Return tag of mask'''
        return self.name

    def getBoardSize(self):
        '''Return the boardSize'''
        return self.boardSize

class Mask_rxcx4(Mask):
    '''Mask that analyzes the states by looking at each row, column, and 4x4 squares.'''

    def __init__(self, name='4x4x4', boardSize=4, maxTile=15):
        '''Init the mask class
        input:
            name: Name of mask used for tag
            boardSize: Size of the board
            maxTile: log2 of max tile that can appear on board'''
        super().__init__(name, boardSize, maxTile)
        # Define some flags to specify the tuple types
        self.row_flag = 0
        self.column_flag = 1
        self.square_flag = 2

    def getNumTuples(self):
        '''Return the number of tuples used to describe each state'''
        return 2*self.boardSize + (self.boardSize-1)**2

    def getMaxTupleNum(self):
        '''Return the largest number corresponding to a tuple'''
        return self.stateToTupleNum(numpy.full((self.boardSize,
                                    self.boardSize), self.maxTile),
                                    self.square_flag,(self.boardSize-1)**2-1)

    def getTupleNums(self, state):
        '''Transforms a state into its tuple number representation
        input:
            state: state to transform
        output: array of tuple nums corresponding to state'''
        tupleNums = numpy.zeros(self.getNumTuples(), dtype=numpy.int)
        index = 0
        # For each tuple type
        for tupleType in range(self.square_flag+1):
            # If tuple type is row or column there will be boardSize tuples
            if tupleType < self.square_flag:
                maxNum = self.boardSize
            # If tuple type is square there will be (self.boardSize - 1)**2
            # tuples
            else:
                maxNum = (self.boardSize - 1)**2
            # For each tuple of this type
            for tupleIndex in range(maxNum):
                # Record this tuple num
                tupleNums[index] = self.stateToTupleNum(state, tupleType, tupleIndex)
                # Increment index
                index += 1
        return tupleNums

    def stateToTupleNum(self, state, tupleType, tupleIndex):
        '''Get specific tupleNum for a state. Accomplished by creating a base
        16 string where each digit represents a tile in the tuple and then
        transfoming to a decimal integer.
        input:
            state: State to transform
            tupleType: Type of tuple given by flag
            tupleIndex: Which tuple of this type we are looking for'''
        # The first digit of hex string is the tuple type
        # The seconde digit is the tuple index
        hexString = ''.join(['{:x}'.format(tupleType), '{:x}'.format(tupleIndex)])
        # The next for digits are the values in the tuple
        if(tupleType == self.row_flag):
            hexString += ''.join(['{:x}'.format(state[tupleIndex][col]) for col
                                  in range(self.boardSize)])
        elif(tupleType == self.column_flag):
            hexString += ''.join(['{:x}'.format(state[row][tupleIndex]) for row
                                  in range(self.boardSize)])
        elif(tupleType == self.square_flag):
            basex = tupleIndex % (self.boardSize - 1)
            basey = int(tupleIndex / (self.boardSize - 1))
            hexString += ''.join(['{:x}'.format(state[basex + i][basey + j])
                                  for i in range(2) for j in range(2)])
        # Return the hexString as as a decimal integer
        return int(hexString, base=16)


text:
#!/usr/bin/env python

from numpy import ones, sqrt

def ispandigital(n):
    s = str(n)
    l = len(s)
    s = s.replace('0', '')  # remove 0
    s = ''.join(sorted(set(s))) # remove duplicates
    if len(s) == l and int(s[-1]) == l:
        return True
    else:
        return False

def test_ispandigital():
    assert ispandigital(12304) == False
    assert ispandigital(1245) == False
    assert ispandigital(12245) == False
    assert ispandigital(123) == True

N = 9999999 + 1
p = ones(N, int)
p[0] = p[1] = 0
for n in range(2, int(sqrt(N))):
    for i in range(2 * n, N, n):
        p[i] = 0

for i in range(N - 1, 0, -1):
    if p[i] == 1:
        if ispandigital(i) == True:
            print i
            break

max_stars_repo_name: autozimu/projecteuler

max_issues_repo_name: autozimu/projecteuler

max_forks_repo_name: autozimu/projecteuler