text:
import numpy as np
import code
import theano
from theano import config
import theano.tensor as tensor
from theano.ifelse import ifelse
from collections import OrderedDict
import time
from imagernn.utils import zipp, initwTh, numpy_floatX, _p, sliceT, basic_lstm_layer, dropout_layer, gumbel_softmax_sample, l2norm
#from theano.tensor.shared_randomstreams import RandomStreams as RandomStreams
from theano.sandbox.rng_mrg import MRG_RandomStreams as RandomStreams

class LSTMGenerator:
  """
  A multimodal long short-term memory (LSTM) generator
  """
# ========================================================================================
  def __init__(self, params):

    image_encoding_size = params.get('image_encoding_size', 128)
    word_encoding_size = params.get('word_encoding_size', 128)

    hidden_size = params.get('hidden_size', 128)
    hidden_depth = params.get('hidden_depth', 1)
    generator = params.get('generator', 'lstm')
    vocabulary_size = params.get('vocabulary_size',-1)
    output_size = params.get('output_size',-1)
    image_feat_size = params.get('image_feat_size',-1)# size of CNN vectors hardcoded here

    aux_inp_size = params.get('aux_inp_size', -1)

    model = OrderedDict()
    # Recurrent weights: take x_t, h_{t-1}, and bias unit
    # and produce the 3 gates and the input to cell signal

    encoder = params.get('feat_encoder', None)
    use_feat_enc = params.get('use_encoder_for',0)

    if not (use_feat_enc & 1):
        model['WIemb'] = initwTh(image_feat_size, word_encoding_size) # image encoder
        model['b_Img'] = np.zeros((word_encoding_size)).astype(config.floatX)

    model['Wemb'] = initwTh(vocabulary_size, word_encoding_size) # word encoder
    model['lstm_W_hid'] = initwTh(hidden_size, 4 * hidden_size)
    model['lstm_W_inp'] = initwTh(word_encoding_size, 4 * hidden_size)

    for i in xrange(1,hidden_depth):
        model['lstm_W_hid_'+str(i)] = initwTh(hidden_size, 4 * hidden_size)
        model['lstm_W_inp_'+str(i)] = initwTh(hidden_size, 4 * hidden_size)

    model['lstm_b'] = np.zeros((4 * hidden_size,)).astype(config.floatX)
    # Decoder weights (e.g. mapping to vocabulary)

    if params.get('class_out_factoring',0) == 0:
        model['Wd'] = initwTh(hidden_size, output_size) # decoder
        model['bd'] = np.zeros((output_size,)).astype(config.floatX)
    else:
        clsinfo = params['ixtoclsinfo']
        self.clsinfo = clsinfo
        clsSizes = clsinfo[:,2] - clsinfo[:,1]
        self.clsSize = np.zeros(params['nClasses'])
        self.clsOffset = np.zeros(params['nClasses'],dtype=np.int32)
        self.clsSize[clsinfo[:,0]] = clsSizes
        self.clsOffset[clsinfo[:,0]] = np.int32(clsinfo[:,1])
        max_cls_size = np.max(clsSizes)
        self.max_cls_size = max_cls_size
        Wd = np.zeros((params['hidden_size'],params['nClasses'], max_cls_size),dtype=config.floatX)
        model['bd'] = np.zeros((1,params['nClasses'], max_cls_size),dtype=config.floatX)
        for cix in clsinfo[:,0]:
            Wd[:,cix,:clsSizes[cix]] = initwTh(params['hidden_size'],clsSizes[cix])
            model['bd'][0,cix,clsSizes[cix]:] = -100
        model['Wd'] = Wd

    update_list = ['lstm_W_hid', 'lstm_W_inp', 'lstm_b', 'Wd', 'bd','Wemb']
    self.regularize = ['lstm_W_hid', 'lstm_W_inp', 'Wd','Wemb']
    if not (use_feat_enc & 1):
        update_list.extend(['WIemb', 'b_Img'])
        self.regularize.extend(['WIemb'])

    if params.get('class_out_factoring',0) == 1:
        model['WdCls'] = initwTh(hidden_size, params['nClasses']) # decoder
        model['bdCls'] = np.zeros((params['nClasses'],)).astype(config.floatX)
        update_list.extend(['WdCls', 'bdCls'])
        self.regularize.extend(['WdCls'])

    for i in xrange(1,hidden_depth):
        update_list.append('lstm_W_hid_'+str(i))
        update_list.append('lstm_W_hid_'+str(i))
        self.regularize.append('lstm_W_inp_'+str(i))
        self.regularize.append('lstm_W_inp_'+str(i))

    if params.get('en_aux_inp',0):
        if params.get('swap_aux',1) == 1:
            if not (use_feat_enc & 2) or params.get('encode_gt_sentences',0):
                model['WIemb_aux'] = initwTh(aux_inp_size, image_encoding_size) # image encoder
                model['b_Img_aux'] = np.zeros((image_encoding_size)).astype(config.floatX)
                update_list.append('WIemb_aux')
                self.regularize.append('WIemb_aux')
                update_list.append('b_Img_aux')
            model['lstm_W_aux'] = initwTh(image_encoding_size, 4 * hidden_size, 0.00005)
        else:
            model['lstm_W_aux'] = initwTh(aux_inp_size, 4 * hidden_size, 0.001)
        update_list.append('lstm_W_aux')
        self.regularize.append('lstm_W_aux')

    if params.get('gen_input_noise',0):
        self.noise_dim = params.get('gen_inp_noise_dim',50)
        model['lstm_W_noise'] = initwTh(self.noise_dim, 4 * hidden_size, 0.001)

    self.model_th = self.init_tparams(model)
    del model
    if params.get('use_gumbel_mse',0):
        self.usegumbel = theano.shared(1)
        self.gumb_temp = theano.shared(numpy_floatX(params.get('gumbel_temp_init',0.5)))
        #self.model_th['gumb_temp'] = self.gumb_temp
        self.softmax_smooth_factor = theano.shared(numpy_floatX(params.get('softmax_smooth_factor',1.0)))
    else:
        self.usegumbel = theano.shared(0)
    self.update_list = update_list

# ========================================================================================
  def init_tparams(self,params):
    tparams = OrderedDict()
    for kk, pp in params.iteritems():
        tparams[kk] = theano.shared(params[kk], name=kk)
    return tparams

# ========================================================================================
 # BUILD LSTM forward propogation model
  def build_model(self, tparams, options, xI=None, xAux = None, attn_nw = None):
    self.trng = RandomStreams(int(time.time()))

    # Used for dropout.
    use_noise = theano.shared(numpy_floatX(0.))

    xW = tensor.matrix('xW', dtype='int64')

    mask = tensor.matrix('mask', dtype=config.floatX)
    n_timesteps = xW.shape[0]
    n_samples = xW.shape[1]

    embW = tparams['Wemb'][xW.flatten()].reshape([n_timesteps,
                                                n_samples,
                                                options['word_encoding_size']])
    if xI == None:
        xI = tensor.matrix('xI', dtype=config.floatX)
        embImg = (tensor.dot(xI, tparams['WIemb']) + tparams['b_Img'])
        xI_is_inp = True
    else:
        embImg = xI
        xI_is_inp = False

    if xAux == None:
        xAux = tensor.matrix('xAux', dtype=config.floatX) if attn_nw == None else tensor.tensor3('xAux', dtype=config.floatX)
        if (options.get('swap_aux',1)) and (attn_nw == None):
           xAuxEmb = tensor.dot(xAux,tparams['WIemb_aux']) + tparams['b_Img_aux']
        else:
           xAuxEmb = xAux
        xA_is_inp = True
    else:
        xA_is_inp = False
        if options.get('encode_gt_sentences',0):
           xAuxEmb = tensor.dot(xAux,tparams['WIemb_aux']) + tparams['b_Img_aux']
        else:
            xAuxEmb = xAux

    embImg = embImg.reshape([1,n_samples,options['image_encoding_size']]);
    emb = tensor.concatenate([embImg, embW], axis=0)

    #This is implementation of input dropout !!
    if options['use_dropout']:
        emb = dropout_layer(emb, use_noise, self.trng, options['drop_prob_encoder'], shp = emb.shape)
        if (options.get('en_aux_inp',0)) and (attn_nw == None):
            xAuxEmb = dropout_layer(xAuxEmb, use_noise, self.trng, options['drop_prob_aux'], shp = xAuxEmb.shape)

    # Implement scehduled sampling!
    if options.get('sched_sampling_mode',None) != None:
        curr_epoch = tensor.scalar(name='curr_epoch',dtype=config.floatX)

        # Assign the probabilies according to the scheduling mode
        if options['sched_sampling_mode'] == 'linear':
            prob = tensor.maximum(options['sslin_min'],options['sched_sampling_const'] - options['sslin_slope'] * curr_epoch)
        elif options['sched_sampling_mode'] == 'exp':
            raise ValueError('ERROR: %s --> This solver type is not yet supported'%(options['sched_sampling_mode']))
        elif options['sched_sampling_mode'] == 'invsig':
            raise ValueError('ERROR: %s --> This solver type is not yet supported'%(options['sched_sampling_mode']))
        else:
            raise ValueError('ERROR: %s --> This scheduling type is unknown'%(options['sched_sampling_mode']))

        # Now to build the mask. We don't want to do this coin toss when
        # feeding in image feature and the start symbol
        sched_mask = self.trng.binomial((n_timesteps - 2, n_samples), p=prob, n=1, dtype='int64')
        sched_mask = tensor.concatenate([sched_mask, tensor.alloc(1, 2, n_samples)],axis=0)
    else:
        sched_mask = []

    #############################################################################################################################
    # This implements core lstm
    rval, updatesLSTM = basic_lstm_layer(tparams, emb[:n_timesteps,:,:], xAuxEmb, use_noise, options,
                                         prefix=options['generator'], sched_prob_mask = sched_mask, attn_nw = attn_nw)
    #############################################################################################################################

    # NOTE1: we are leaving out the first prediction, which was made for the image and is meaningless.
    if options['use_dropout']:
        # XXX : Size given to dropout is missing one dimension. This keeps the dropped units consistent across time!?.
        # ###   Is this a good bug ?
        p = dropout_layer(sliceT(rval[0][1:,:,:],options.get('hidden_depth',1),options['hidden_size']), use_noise, self.trng,
            options['drop_prob_decoder'], (n_samples,options['hidden_size']))
    else:
        p = sliceT(rval[0][1:,:,:],options.get('hidden_depth',1),options['hidden_size'])

    if options.get('class_out_factoring',0) == 1:
        if options.get('cls_diff_layer',0) == 1:
            pC_inp = dropout_layer(sliceT(rval[0][1:,:,:],options.get('hidden_depth',1)-2,options['hidden_size']), use_noise, self.trng,
            options['drop_prob_decoder'], (n_samples,options['hidden_size']))
        else:
            pC_inp = p

    n_out_samps = (n_timesteps-1) * n_samples
    if options.get('class_out_factoring',0) == 0:
        pW = (tensor.dot(p,tparams['Wd']) + tparams['bd']).reshape([n_out_samps,options['output_size']])
        if options.get('use_gumbel_mse',0) == 0:
            pWSft = tensor.nnet.softmax(pW)
        else:
            w_out = ifelse(self.usegumbel, gumbel_softmax_sample(self.trng, pW, self.gumb_temp,
                            hard=options.get('use_gumbel_hard',False)), tensor.nnet.softmax(pW))
            # This is not exactly right, but just testing
            pWSft = w_out

        totProb = pWSft[tensor.arange(n_out_samps), xW[1:,:].flatten()]
        out_list = [pWSft, totProb, pW]
    else:
        ixtoclsinfo_t = tensor.as_tensor_variable(self.clsinfo)
        xC = ixtoclsinfo_t[xW[1:,:].flatten(),0]
        if options.get('cls_zmean',1):
            pW = ((tparams['Wd'][:,xC,:].T*((p.reshape([1,n_out_samps,options['hidden_size']])-tparams['WdCls'][:,xC].T))).sum(axis=-1).T
                 + tparams['bd'][:,xC,:])
        else:
            pW = ((tparams['Wd'][:,xC,:].T*(p.reshape([1,n_out_samps,options['hidden_size']]))).sum(axis=-1).T
                 + tparams['bd'][:,xC,:])
        pWSft   = tensor.nnet.softmax(pW[0,:,:])

        pC    = (tensor.dot(pC_inp,tparams['WdCls']) + tparams['bdCls']).reshape([n_out_samps,options['nClasses']])
        pCSft = tensor.nnet.softmax(pC)

        totProb = pWSft[tensor.arange(n_out_samps), ixtoclsinfo_t[xW[1:,:].flatten(),3]] * \
                  pCSft[tensor.arange(n_out_samps), xC]
        out_list = [pWSft, pCSft, totProb, p]

    tot_cost = -(tensor.log(totProb + 1e-10) * mask[1:,:].flatten()).sum()
    tot_pplx = -(tensor.log2(totProb + 1e-10) * mask[1:,:].flatten()).sum()
    cost = [tot_cost/tensor.cast(n_samples,dtype = config.floatX), tot_pplx]

    inp_list = [xW, mask]
    if xI_is_inp:
        inp_list.append(xI)

    if options.get('en_aux_inp',0) and xA_is_inp:
        inp_list.append(xAux)

    if options.get('sched_sampling_mode',None) != None:
        inp_list.append(curr_epoch)

    f_pred_prob = theano.function([xW, xI, xAux], out_list, name='f_pred_prob', updates=updatesLSTM)

    return use_noise, inp_list, f_pred_prob, cost, out_list , updatesLSTM

# ========================================================================================
# Predictor Related Stuff!!

  def prepPredictor(self, model_npy=None, checkpoint_params=None, beam_size=5, xI=None, xAux = None, inp_list_prev=[], per_word_logweight = None):
    if model_npy != None:
        if type(model_npy[model_npy.keys()[0]]) == np.ndarray:
            zipp(model_npy, self.model_th)
        else:
            self.model_th = model_npy

    #theano.config.exception_verbosity = 'high'
    self.beam_size = beam_size

	# Now we build a predictor model
    if checkpoint_params.get('advers_gen',0) == 1:
        checkpoint_params['n_gen_samples'] = beam_size
    (inp_list_gen, predLogProb, predIdx, predCand, wOut_emb, updates, seq_lengths) = self.build_prediction_model(self.model_th, checkpoint_params, xI, xAux, per_word_logweight = per_word_logweight)
    self.f_pred_th = theano.function(inp_list_prev + inp_list_gen, [predLogProb, predIdx, predCand], name='f_pred')

	# Now we build a training model which evaluates cost. This is for the evaluation part in the end
    if checkpoint_params.get('advers_gen',0) == 0:
        (self.use_dropout, inp_list_gen2,
         f_pred_prob, cost, predTh, updatesLSTM) = self.build_model(self.model_th, checkpoint_params, xI, xAux)
        self.f_eval= theano.function(inp_list_prev+inp_list_gen2, cost, name='f_eval')

# ========================================================================================
  def predict(self, batch, checkpoint_params, ext_inp = []):

    inp_list = ext_inp
    if not checkpoint_params.get('use_encoder_for',0)&1:
        inp_list.extend([batch[0]['image']['feat'].reshape(1,checkpoint_params['image_feat_size']).astype(config.floatX)])

    if not checkpoint_params.get('use_encoder_for',0)&2:
        if checkpoint_params.get('en_aux_inp',0):
            inp_list.append(batch[0]['image']['aux_inp'].reshape(1,checkpoint_params['aux_inp_size']).astype(config.floatX))

    Ax = self.f_pred_th(*inp_list)

    # Backtracking to decode the correct sequence of candidates
    Ys = []
    for i in xrange(self.beam_size):
        candI = []
        curr_cand = Ax[2][-1][i]
        for j in reversed(xrange(Ax[1].shape[0]-1)):
            candI.insert(0,Ax[1][j][curr_cand])
            curr_cand = Ax[2][j][curr_cand]

        Ys.append([Ax[0][i], candI])
    return [Ys], Ax

  def build_prediction_model(self, tparams, options, xI=None, xAux = None, per_word_logweight = None):
    #Initialize random streams for other to use.
    self.trng = RandomStreams(int(time.time()))

    if xI == None:
        xI = tensor.matrix('xI', dtype=config.floatX)
        embImg = (tensor.dot(xI, tparams['WIemb']) + tparams['b_Img'])
        xI_is_inp = True
    else:
        xI_is_inp = False
        embImg = xI
    if xAux == None and options.get('en_aux_inp',0):
        xAux = tensor.matrix('xAux', dtype=config.floatX)
        xA_is_inp = True
        if options.get('swap_aux',1):
           xAuxEmb = tensor.dot(xAux,tparams['WIemb_aux']) + tparams['b_Img_aux']
        else:
           xAuxEmb = xAux
    else:
        xA_is_inp = False
        if options.get('encode_gt_sentences',0):
           xAuxEmb = tensor.dot(xAux,tparams['WIemb_aux']) + tparams['b_Img_aux']
        else:
           xAuxEmb = xAux

    if options.get('advers_gen',0) == 1:
        accLogProb, Idx, wOut_emb, updates, seq_lengths = self.lstm_advers_gen_layer(tparams, embImg, xAuxEmb, options, prefix=options['generator'])
        Cand = tensor.tile(tensor.arange(Idx.shape[1]),[Idx.shape[0],1])
    else:
        accLogProb, Idx, Cand, wOut_emb, updates = self.lstm_predict_layer(tparams, embImg, xAuxEmb, options, self.beam_size, prefix=options['generator'],
                per_word_logweight=per_word_logweight)
        seq_lengths = []

    inp_list = []
    if xI_is_inp:
        inp_list.append(xI)
    if options.get('en_aux_inp',0) and xA_is_inp:
        inp_list.append(xAux)

    return inp_list, accLogProb, Idx, Cand, wOut_emb, updates, seq_lengths

# ========================================================================================
  # LSTM LAYER in Prediction mode. Here we don't provide the word sequences, just the image feature vector
  # The network starts first with forward propogatin the image feature vector. Then we pass the start word feature
  # i.e zeroth word vector. From then the network output word (i.e ML word) is fed as the input to the next time step.
  # In beam_size > 1 we could repeat a time step multiple times, once for each beam!!.

  def lstm_predict_layer(self, tparams, Xi, aux_input, options, beam_size, prefix='lstm', per_word_logweight = None):

    nMaxsteps = options.get('maxlen',30)

    if nMaxsteps is None:
        nMaxsteps = 30
    n_samples = 1
    h_depth = options.get('hidden_depth',1)
    h_sz = options['hidden_size']

    # ----------------------  STEP FUNCTION  ---------------------- #
    def _stepP(x_, h_, c_, lP_, dV_, xAux):
        preact = tensor.dot(sliceT(h_, 0, h_sz), tparams[_p(prefix, 'W_hid')])
        preact += (tensor.dot(x_, tparams[_p(prefix, 'W_inp')]) +
                   tparams[_p(prefix, 'b')])
        if options.get('en_aux_inp',0):
            preact += tensor.dot(xAux,tparams[_p(prefix,'W_aux')])

        hL = [[]]*h_depth
        cL = [[]]*h_depth
        outp = [[]]*h_depth
        for di in xrange(h_depth):
            i = tensor.nnet.sigmoid(sliceT(preact, 0, h_sz))
            f = tensor.nnet.sigmoid(sliceT(preact, 1, h_sz))
            o = tensor.nnet.sigmoid(sliceT(preact, 2, h_sz))
            cL[di] = tensor.tanh(sliceT(preact, 3, h_sz))
            cL[di] = f * sliceT(c_, di, h_sz) + i * cL[di]
            hL[di] = o * tensor.tanh(cL[di])
            outp[di] = hL[di]
            if options.get('en_residual_conn',1):
                if (di > 0):
                  outp[di] += outp[di-1]
                  print "Connecting residual at %d"%(di)
            if di < (h_depth - 1):
                preact = tensor.dot(sliceT(h_, di+1, h_sz), tparams[_p(prefix, ('W_hid_' + str(di+1)))]) + \
                        tensor.dot(outp[di], tparams[_p(prefix, ('W_inp_' + str(di+1)))])

        c = tensor.concatenate(cL,axis=1)
        h = tensor.concatenate(hL,axis=1)

        if options.get('class_out_factoring',0) == 1:
            if options.get('cls_diff_layer',0) == 1:
                pC    = tensor.dot(hL[-2],tparams['WdCls']) + tparams['bdCls']
            else:
                pC    = tensor.dot(outp[-1],tparams['WdCls']) + tparams['bdCls']

            pCSft = tensor.nnet.softmax(pC)
            xCIdx =  tensor.argmax(pCSft,axis=-1)
            #pW = tensor.dot(outp[-1],tparams['Wd'][:,xCIdx,:]) + tparams['bd'][:,xCIdx,:]
            #smooth_factor = tensor.as_tensor_variable(numpy_floatX(options.get('softmax_smooth_factor',1.0)), name='sm_f')
            #pWSft = tensor.nnet.softmax(pW*smooth_factor)
            #lProb = tensor.log(pWSft + 1e-20) + tensor.log(pCSft[0,xCIdx] + 1e-20)
            #########################################################
            # pW is now of size (beam_size, n_classes, class_size)
            if options.get('cls_zmean',0):
                pW = tensor.dot((outp[-1]-tparams['WdCls'][:,xCIdx].T),tparams['Wd'].swapaxes(0,1)) + tparams['bd'][0,:,:]
            else:
                pW = tensor.dot((outp[-1]),tparams['Wd'].swapaxes(0,1)) + tparams['bd'][0,:,:]
            #smooth_factor = tensor.as_tensor_variable(numpy_floatX(options.get('softmax_smooth_factor',1.0)), name='sm_f')
            pWSft = tensor.nnet.softmax(pW.reshape([pW.shape[0]*pW.shape[1],
                    pW.shape[2]])).reshape([pW.shape[0], pW.shape[1]*pW.shape[2]])
            ixtoclsinfo_t = tensor.as_tensor_variable(self.clsinfo)
            lProb = tensor.log(pWSft[:,ixtoclsinfo_t[:,0]*tparams['Wd'].shape[2]+ixtoclsinfo_t[:,3]] + 1e-20) + \
                    tensor.log(pCSft[0,ixtoclsinfo_t[:,0]] + 1e-20)
        else:
            p = tensor.dot(outp[-1],tparams['Wd']) + tparams['bd']
            smooth_factor = tensor.as_tensor_variable(numpy_floatX(options.get('softmax_smooth_factor',1.0)), name='sm_f')
            p = tensor.nnet.softmax(p*smooth_factor)
            lProb = tensor.log(p + 1e-20)
            if per_word_logweight is not None:
                log_w = theano.shared(per_word_logweight)#, dtype= theano.config.floatX)
                lProb = log_w+lProb

        if beam_size > 1:
            def _FindB_best(lPLcl, lPprev, dVLcl):
                srtLcl = tensor.argsort(-lPLcl)
                srtLcl = srtLcl[:beam_size]
                deltaVec = tensor.fill( lPLcl[srtLcl], numpy_floatX(-10000.))
                deltaVec = tensor.set_subtensor(deltaVec[0], lPprev)
                lProbBest = ifelse(tensor.eq( dVLcl, tensor.zeros_like(dVLcl)), lPLcl[srtLcl] + lPprev, deltaVec)
                xWIdxBest = ifelse(tensor.eq( dVLcl, tensor.zeros_like(dVLcl)), srtLcl, tensor.zeros_like(srtLcl))
                return lProbBest, xWIdxBest

            rvalLcl, updatesLcl = theano.scan(_FindB_best, sequences = [lProb, lP_, dV_], name=_p(prefix, 'FindBest'), n_steps=x_.shape[0])
            xWIdxBest = rvalLcl[1]
            lProbBest = rvalLcl[0]

            xWIdxBest = xWIdxBest.flatten()
            lProb = lProbBest.flatten()
            # Now sort and find the best among these best extensions for the current beams
            srtIdx = tensor.argsort(-lProb)
            srtIdx = srtIdx[:beam_size]
            xCandIdx = srtIdx // beam_size # Floor division
            h = h.take(xCandIdx.flatten(),axis=0)
            c = c.take(xCandIdx.flatten(),axis=0)
            xWlogProb = lProb[srtIdx]
            xWIdx = xWIdxBest[srtIdx]
            if options.get('class_out_factoring',0) == 1:
                clsoffset = tensor.as_tensor_variable(self.clsOffset)
        else:
            xCandIdx = tensor.as_tensor_variable([0])
            lProb = lProb.flatten()
            xWIdx =  tensor.argmax(lProb,keepdims=True)
            xWlogProb = lProb[xWIdx] + lP_
            #if options.get('class_out_factoring',0) == 1:
            #    clsoffset = tensor.as_tensor_variable(self.clsOffset)
            #    xWIdx += clsoffset[xCIdx]
            h = h.take(xCandIdx.flatten(),axis=0)
            c = c.take(xCandIdx.flatten(),axis=0)

        if options.get('softmax_propogate',0) == 0:
            xW = tparams['Wemb'][xWIdx.flatten()]
        else:
            xW = p.dot(tparams['Wemb'])
        doneVec = tensor.eq(xWIdx,tensor.zeros_like(xWIdx))

        return [xW, h, c, xWlogProb, doneVec, xWIdx, xCandIdx], theano.scan_module.until(doneVec.all())
    # ------------------- END of STEP FUNCTION  -------------------- #

    if options.get('en_aux_inp',0) == 0:
       aux_input = []

    h = tensor.alloc(numpy_floatX(0.),beam_size,h_sz*h_depth)
    c = tensor.alloc(numpy_floatX(0.),beam_size,h_sz*h_depth)

    lP = tensor.alloc(numpy_floatX(0.), beam_size);
    dV = tensor.alloc(np.int8(0.), beam_size);

    # Propogate the image feature vector
    [xW, h, c, _, _, _, _], _ = _stepP(Xi, h[:1,:], c[:1,:], lP, dV,aux_input)

    xWStart = tparams['Wemb'][[0]]
    [xW, h, c, lP, dV, idx0, cand0], _ = _stepP(xWStart, h[:1,:], c[:1,:], lP, dV, aux_input)

    if options.get('en_aux_inp',0) == 1:
        aux_input = tensor.extra_ops.repeat(aux_input,beam_size,axis=0)

    # Now lets do the loop.
    rval, updates = theano.scan(_stepP, outputs_info=[xW, h, c, lP, dV, None, None], non_sequences = [aux_input], name=_p(prefix, 'predict_layers'), n_steps=nMaxsteps)

    return rval[3][-1], tensor.concatenate([idx0.reshape([1,beam_size]), rval[5]],axis=0), tensor.concatenate([cand0.reshape([1,beam_size]), rval[6]],axis=0), tensor.concatenate([tensor.shape_padleft(xW,n_ones=1),rval[0]],axis=0), updates

#================================================================================================================
  def lstm_advers_gen_layer(self, tparams, xI, xAux, options, prefix='lstm'):
    nBatchSamps = xI.shape[0]
    nMaxsteps = options.get('maxlen',15)
    if nMaxsteps is None:
        nMaxsteps = 30
    n_samp = options.get('n_gen_samples',1)

    h_depth = options.get('hidden_depth',1)
    h_sz = options['hidden_size']

    # ----------------------  STEP FUNCTION  ---------------------- #
    def _stepP(U, xW_, h_, c_, lP_, dV_, xAux, xNoise):
        preact = tensor.dot(sliceT(h_, 0, h_sz), tparams[_p(prefix, 'W_hid')])
        preact += (tensor.dot(xW_, tparams[_p(prefix, 'W_inp')]) +
                   tparams[_p(prefix, 'b')])
        preact += xAux
        if options.get('gen_input_noise',0):
            preact += xNoise

        hL = [[]]*h_depth
        cL = [[]]*h_depth
        outp = [[]]*h_depth
        for di in xrange(h_depth):
            i = tensor.nnet.sigmoid(sliceT(preact, 0, h_sz))
            f = tensor.nnet.sigmoid(sliceT(preact, 1, h_sz))
            o = tensor.nnet.sigmoid(sliceT(preact, 2, h_sz))
            cL[di] = tensor.tanh(sliceT(preact, 3, h_sz))
            cL[di] = f * sliceT(c_, di, h_sz) + i * cL[di]
            hL[di] = o * tensor.tanh(cL[di])
            outp[di] = hL[di]
            if options.get('en_residual_conn',1):
                if (di > 0):
                  outp[di] += outp[di-1]
                  print "Connecting residual at %d"%(di)
            if di < (h_depth - 1):
                preact = tensor.dot(sliceT(h_, di+1, h_sz), tparams[_p(prefix, ('W_hid_' + str(di+1)))]) + \
                        tensor.dot(outp[di], tparams[_p(prefix, ('W_inp_' + str(di+1)))])

        c = tensor.concatenate(cL,axis=1)
        h = tensor.concatenate(hL,axis=1)

        logits = tensor.dot(outp[-1],tparams['Wd']) + tparams['bd']
        #p = tensor.dot(outp[-1],l2norm(tparams['Wd'],axis=0))# + tparams['bd']
        if options.get('use_gumbel_mse',0) == 0 or options.get('greedy',0):
            p = tensor.nnet.softmax(logits)
        else:
            p = gumbel_softmax_sample(self.trng, logits*self.softmax_smooth_factor, self.gumb_temp, U, options.get('use_gumbel_hard',False))

        if options.get('computelogprob', 0):
            lProb = tensor.log(tensor.nnet.softmax(logits*self.softmax_smooth_factor)+1e-20)
        else:
            lProb = logits

        # Idx of the correct word should come from the
        xWIdx =  ~dV_ * tensor.argmax(p,axis = -1)

        xWlogProb = ~dV_ * lProb[tensor.arange(nBatchSamps*n_samp),xWIdx] + lP_
        #xW = tparams['Wemb'][xWIdx.flatten()]
        if options.get('use_gumbel_hard',0) and options.get('use_gumbel_mse', 0) and not options.get('greedy', 0):
            xW = p.dot(tparams['Wemb'])
        else:
            xW = theano.gradient.disconnected_grad(tparams['Wemb'][xWIdx.flatten()].reshape([xWIdx.shape[0],-1]))

        doneVec = tensor.eq(xWIdx,tensor.zeros_like(xWIdx))

        return [xW, h, c, xWlogProb, doneVec, xWIdx, p], theano.scan_module.until(doneVec.all())
    # ------------------- END of STEP FUNCTION  -------------------- #

    if options.get('use_gumbel_mse',0) == 0:
        U = self.trng.uniform((nMaxsteps, 1),low=0.,
                              high=1.,dtype=theano.config.floatX)
    else:
        U = self.trng.uniform((nMaxsteps+1, nBatchSamps * n_samp, options['vocabulary_size']),low=0.,
                              high=1.,dtype=theano.config.floatX)

    xI= tensor.extra_ops.repeat(xI, n_samp,axis=0)
    xAux = tensor.extra_ops.repeat(tensor.dot(xAux,
                tparams[_p(prefix,'W_aux')]), n_samp,axis=0)

    if options.get('gen_input_noise',0):
        xNoise = tensor.dot(self.trng.normal([nBatchSamps * n_samp, self.noise_dim]), tparams[_p(prefix,'W_noise')])
    else:
        xNoise = []

    if options.get('gen_use_rand_init',0) and not options.get('gen_input_noise',0):
        h = tensor.unbroadcast(self.trng.uniform([nBatchSamps * n_samp, h_sz*h_depth],low=-0.1, high= 0.1),0,1)
        c = tensor.unbroadcast(self.trng.uniform([nBatchSamps * n_samp,h_sz*h_depth],low=-0.1, high= 0.1),0,1)
    else:
        h = tensor.zeros([nBatchSamps * n_samp, h_sz*h_depth])
        c = tensor.zeros([nBatchSamps * n_samp,h_sz*h_depth])

    lP = tensor.alloc(numpy_floatX(0.), nBatchSamps * n_samp);
    dV = tensor.alloc(np.bool_(0.), nBatchSamps * n_samp);

    # Propogate the image feature vector
    [_, h, c, _, _, _, _], _ = _stepP(U[0,:], xI, h, c, lP, dV, xAux, xNoise)

    xWStart = tensor.unbroadcast(tensor.tile(tparams['Wemb'][[0]],[nBatchSamps * n_samp,1]) ,0,1)

    # Now lets do the loop.
    rval, updates = theano.scan(_stepP, sequences = [U[1:,:]], outputs_info=[xWStart, h, c, lP, dV, None, None],
                                        non_sequences = [xAux, xNoise], name=_p(prefix, 'adv_predict_layers'),
                                        n_steps=nMaxsteps)

    seq_lengths = theano.gradient.disconnected_grad(tensor.argmax(tensor.concatenate([rval[4][:-1,:],
                    tensor.ones((1,nBatchSamps * n_samp))], axis=0), axis=0)+1)

    return rval[3][-1], rval[5], rval[6], updates, seq_lengths

# ========================================================================================
  def build_eval_other_sent(self, tparams, options,model_npy):

    zipp(model_npy, self.model_th)

    # Used for dropout.
    use_noise = theano.shared(numpy_floatX(0.))

    xW = tensor.matrix('xW', dtype='int64')
    mask = tensor.matrix('mask', dtype=config.floatX)
    n_timesteps = xW.shape[0]
    n_samples = xW.shape[1]
    n_out_samps = (n_timesteps-1) * n_samples

    embW = tparams['Wemb'][xW.flatten()].reshape([n_timesteps,
                                                n_samples,
                                                options['word_encoding_size']])
    xI = tensor.matrix('xI', dtype=config.floatX)
    xAux = tensor.matrix('xAux', dtype=config.floatX)

    if options.get('swap_aux',0):
       xAuxEmb = tensor.dot(xAux,tparams['WIemb_aux']) + tparams['b_Img_aux']
    else:
       xAuxEmb = xAux

    embImg = (tensor.dot(xI, tparams['WIemb']) + tparams['b_Img']).reshape([1,n_samples,options['image_encoding_size']]);
    emb = tensor.concatenate([embImg, embW], axis=0)

    rval, updatesLSTM = basic_lstm_layer(tparams, emb[:n_timesteps,:,:], xAuxEmb, use_noise, options, prefix=options['generator'])
    p = sliceT(rval[0][1:,:,:],options.get('hidden_depth',1),options['hidden_size'])

    if options.get('class_out_factoring',0) == 0:
        pW = (tensor.dot(p,tparams['Wd']) + tparams['bd']).reshape([n_out_samps,options['output_size']])
        pWSft = tensor.nnet.softmax(pW)
        totProb = pWSft[tensor.arange(n_out_samps), xW[1:,:].flatten()]
        out_list = [pWSft, totProb, p]
    else:
        ixtoclsinfo_t = tensor.as_tensor_variable(self.clsinfo)
        xC = ixtoclsinfo_t[xW[1:,:].flatten(),0]
        pW = ((tparams['Wd'][:,xC,:].T*((p.reshape([1,n_out_samps,options['hidden_size']])-tparams['WdCls'][:,xC].T))).sum(axis=-1).T
             + tparams['bd'][:,xC,:])
        pWSft   = tensor.nnet.softmax(pW[0,:,:])
        pC    = (tensor.dot(p,tparams['WdCls']) + tparams['bdCls']).reshape([n_out_samps,options['nClasses']])
        pCSft = tensor.nnet.softmax(pC)

        totProb = pWSft[tensor.arange(n_out_samps), ixtoclsinfo_t[xW[1:,:].flatten(),3]] * \
                  pCSft[tensor.arange(n_out_samps), xC]

    tot_cost = -(tensor.log(totProb + 1e-10) * mask[1:,:].flatten()).reshape([n_timesteps-1,n_samples])
    cost = tot_cost.sum(axis=0)

    inp_list = [xW, mask, xI]

    if options.get('en_aux_inp',0):
        inp_list.append(xAux)

    self.f_pred_prob_other = theano.function([xW,xI, xAux], pWSft, name='f_pred_prob', updates=updatesLSTM)
    #f_pred = theano.function([xW, mask], pred.argmax(axis=1), name='f_pred')

    #cost = -tensor.log(pred[tensor.arange(n_timesteps),tensor.arange(n_samples), xW] + 1e-8).mean()

    self.f_eval_other = theano.function(inp_list, cost, name='f_eval')

    return use_noise, inp_list, self.f_pred_prob_other, cost, pW, updatesLSTM

# =================================== MULTI Model Ensemble Predictor related ==========================

max_stars_repo_name: rakshithShetty/captionGAN

reason: import numpy,import theano,from theano


text:
"""Synthesize coherent anti-stokes Raman spectra for common gas species."""
import numpy as np
from .line_strength import LineStrength, linewidth_isolated
from ._constants import univ_const, chi_const
from .convol_fcn import gaussian_line, lorentz_line
from .utils import comp_normalize, eq_comp

class CarsSpectrum():
    """Basic class for synthesizing CARS spectrum."""

    def __init__(self, pressure=1, init_comp=None, chi_set="SET 1", **kwargs):
        """Input sample volume conditions.

        Parameters
        ----------
        pressure : float, optional
            Pressure in bars, by default 1.
        init_comp : dict, optional
            Initial composition in the measurement volume (prior to chemical
            reaction). If not given standard air composition is assumed as:

            .. code-block:: python

                {
                'N2': 0.79,
                'Ar': 0.0,
                'CO2': 0,
                'CO': 0,
                'H2': 0,
                'O2': 0.21,
                'H2O': 0
                }
        chi_set : str, optional
            Choose from the available set of non-resonant susceptibilities, by
            default 'SET 1':

            * 'SET 1': based on CARSFT :cite:`Palmer:89`.
            * 'SET 2': based on NRC :cite:`Parameswaran:88`.
            * 'SET 3': based on Eckbreth :cite:`Eckbreth:96`.

        Other Parameters
        ----------------
        species : str, optional
            Specify the molecule, by default 'N2'. Currently only supports
            'N2'.
        custom_dict : dict, optional
            Specify custom molecular constants, by default None. This can be
            used to modify default molecular constants and/or add custom
            species. The dictionary should contain all the necessary keys as in
            the default:

            .. code-block:: python

                ['we', 'wx', 'wy', 'wz', 'Be', 'De', 'alpha_e', 'beta_e',
                'gamma_e', 'H0', 'He', 'mu', 'MW', 'Const_Raman', 'G/A']
        """
        self.pressure = pressure

        if init_comp is None:
            # Assume air
            self.init_comp = {
                'N2': 0.79,
                'Ar': 0.0,
                'CO2': 0,
                'CO': 0,
                'H2': 0,
                'O2': 0.21,
                'H2O': 0,
                }
        else:
            # normalize the mole fraction to 1
            self.init_comp = comp_normalize(init_comp)

        # initiate line strength factors
        self.ls_factors = LineStrength(**kwargs)

        # load Universal constants
        self.Const_N = univ_const("Const_N")  # *1e22 [molecules/cm^3]
        self.C = univ_const("c")  # speed of light

        # load third order nonresonant susceptibilities
        self.chi_nrs = chi_const(chi_set)["SPECIES"]
        self.chi_nrs_T0 = chi_const(chi_set)["T0"]
        self.chi_nrs_P0 = chi_const(chi_set)["P0"]

        # load species-specific constants
        # alpha'/(2*pi*c*v0*m_reduced)^0.5 [*1e-11]
        self.Const_Raman = self.ls_factors.mc_dict['Const_Raman']
        # gamma'/alpha'
        self.pol_ratio = self.ls_factors.mc_dict['G/A']

    def num_dens(self, temperature):
        r"""Calculate number density in the probe volume.

        Parameters
        ----------
        temperature : float
            Temperature in [K].

        Returns
        -------
        float
            Number density in the probe volume in
            [:math:`\mathrm{molecules}/\mathrm{cm}^3`].
        """
        return self.pressure/temperature*self.Const_N

    def chi_nrs_est(self, temperature, local_comp=None):
        """Calculate effective local nonresonant susceptibility.

        Parameters
        ----------
        temperature : float
            Temperature in [K].
        local_comp : dict, optional
            A dictionary of local gas composition (percentile) for the purpose
            of calculating nonresonant background and collisional narrowing, by
            default None. If not given, initial composition will be used.
        """
        if local_comp is None:
            local_comp = self.init_comp
        local_comp = comp_normalize(local_comp)

        # Calculate the effective nonresonant susceptibility based on given
        # composition or air
        _common_keys = [_key for _key in self.chi_nrs if _key in local_comp]
        chi_nrs_eff = np.sum([self.chi_nrs[_key]*local_comp[_key]
                              for _key in _common_keys])

        # Convert to actual pressure and temperature condition
        chi_nrs_eff = chi_nrs_eff*(
            self.pressure/self.chi_nrs_P0*self.chi_nrs_T0/temperature)

        return chi_nrs_eff

    def trans_amp(self, v, j, branch=0):
        """Transition amplitude d or polarizability matrix element.

        It's related to Raman differential scattering cross-section and number
        density. Note that it is molecule based (divided number density N).

        Parameters
        ----------
        v : int
            Vibrational quantum number.
        j : int
            Rotational quantum number.
        branch : int, optional
            Q, S or O-branch with a value of 0, 2 or -2, by default 0.

        Returns
        -------
        d_squared : float
            Transition amplitude squared based on specified v, j and branch.
        """
        # Calculate the line strength factors using the LineStrength class
        pt_coeff, cd_coeff = self.ls_factors.int_corr(j, branch)

        # Calculate the branch-dependent d^2
        if branch in (2, -2):
            d_squared = (self.Const_Raman**2*(4/45)*pt_coeff
                         * self.pol_ratio**2*(v+1)*cd_coeff)
        elif branch == 0:
            d_squared = self.Const_Raman**2*(
                1 + (4/45)*pt_coeff*self.pol_ratio**2)*(v+1)*cd_coeff
        else:
            raise ValueError("Branch can only be 0, 2 or -2")

        return d_squared

    def relax_rate(self, temperature, j_i, j_j, fit_param=None):
        """Calculate relaxation rates.

        This is carried out using the modified exponential gap law
        (MEG) :cite:`Rahn:86`.

        .. note::
            Only considers Q-branch. O- and S-branches are taken as the same as
            Q-branch. It is also independent of vibrational state. Only valid
            for N2 at the moment.

        .. attention::

            Due to nuclear-spin selection rules for N2, gamma_ji is zero when
            (j_j-j_i) is odd. This would need to be adjusted for other types of
            molecules.

        Parameters
        ----------
        temperature : float
            Temperature in [K].
        j_i : int
            Rotational quantum number of the initial state.
        j_j : int
            Rotational quantum number of the final state.
        fit_param : list
            The four fitting parameters in a list: alpha, beta, sigma, m,
            default None. Needs to be provided.

        Returns
        -------
        gamma_ji, gamma_ij : float
            Upward and downward relaxation rate.
        """
        # Energies of the rovibrational states involved.
        Ej_i = self.ls_factors.term_values(0, j_i, 'Fv')
        del_E = self.ls_factors.term_values(0, j_j, 'Fv') - Ej_i

        # For a transition from j_i to j_j with i<j:
        if abs(j_i-j_j) % 2 != 0:
            # For N2 gamma_ji is zero when del_j is odd due to nuclear-spin
            # selection rules.
            gamma_ji = 0
            gamma_ij = 0
        else:
            alpha, beta, sigma, m = fit_param

            _term_1 = (1-np.exp(-m))/(
                1-np.exp(-m*temperature/295))*(295/temperature)**0.5
            _term_2 = ((1+1.5*1.44*Ej_i/temperature/sigma)/(
                1+1.5*1.44*Ej_i/temperature))**2

            gamma_ji = alpha*self.pressure/1.01325*_term_1*_term_2*np.exp(
                -beta*del_E*1.44/temperature)
            gamma_ij = gamma_ji*(2*j_i+1)/(2*j_j+1)*np.exp(
                del_E*1.44/temperature)

        return gamma_ji, gamma_ij

    def relax_mat(self, temperature, js=70, mode='MEG'):
        """Calculate relaxation rate matrix.

        The matrix elements are calculated using different exponential gap
        laws. The calculations are conducted for a fixed number of j levels,
        independent of v or branch.

        .. note::

            The (absolute) diagonal values are the HWHM of Raman linewidth.
            They are cross checked against NRC report :cite:`Parameswaran:88`
            on D4-D9 at 1 atm and 1550 K.

        Parameters
        ----------
        temperature : float
            Temperature in [K].
        js : int, optional
            Total number of rotational levels to be considered, by default 70.
        mode : str, optional
            By default 'MEG' is used. Possible options are:

            * 'MEG': Modified exponential gap law considering only N2-N2
              collisions :cite:`Rahn:86`. The fitting parameters are taken from
              :cite:`Woyde:90`.
            * 'EMEG': Extended MEG which weighs the contributions of N2, CO2
              and H2O :cite:`Woyde:90` (to be implemented).
            * 'XMEG': Extended MEG plus vibrational dephasing :cite:`Porter:90`
              (to be implemented).

        Returns
        -------
        2d array with the shape of (js, js)
            Relaxation matrix with the shape of js x js.
        """
        if mode != 'MEG':
            raise ValueError('Only MEG available at the moment!')
        # Define fitting parameters for N2-N2, N2-CO2, N2-H20
        #               alpha     beta, sigma, m
        #             [cm^-1/atm]
        fit_param_N2 = [0.0231, 1.67, 1.21, 0.1487]

        # Construct the matrix
        gamma_mat = np.zeros([js, js])
        for _i in range(js):
            for _j in range(_i+1, js):  # i.e., _i<_j
                gamma_mat[_j, _i], gamma_mat[_i, _j] = self.relax_rate(
                    temperature, _i, _j, fit_param_N2)

        # Calculate the diagonal
        for _i in range(js):
            gamma_mat[_i, _i] = -np.sum(gamma_mat[:, _i])

        return gamma_mat

    def peak_check(self, x_mol, temperature, v, j, branch=0, Gamma_j=None):
        r"""Calculate the theoretical peak intensity.

        This is done for the specified transition, assuming isolated lines.

        .. note::

            The purpose of this method is to check against the values tabulated
            in the NRC report :cite:`Parameswaran:88` on D4-D9. The agreement
            is excellent if x_mol is set to 1 and if the Gamma_j from NRC is
            used.

        Parameters
        ----------
        x_mol : float
            Mole fraction of probed molecule within [0, 1].
        temperature : float
            Temperature in the probe volume.
        v : int
            Vibrational quantum number.
        j : int
            Rotational quantum number.
        branch : int, optional
            Q, S or O-branch with a value of 0, 2 or -2, by default 0.
        Gamma_j : float, optional
            Raman linewidth, by default None. When not given, the value will be
            calculated from
            :mod:`carspy.line_strength.linewidth_isolated`.

        Returns
        -------
        float
            Peak transition amplitude based on specified v, j and branch.
        """
        # Species number density
        N = x_mol*self.num_dens(temperature)

        # Raman linewidth (depreciated, use MEG instead)
        if Gamma_j is None:
            Gamma_j = linewidth_isolated(self.pressure, temperature, j)

        # Calculate intensity
        del_pop = self.ls_factors.pop_factor(temperature, v, j, branch)
        d_squared = self.trans_amp(v, j, branch)

        return N*del_pop*d_squared/Gamma_j/(2*np.pi*self.C)

    def chi_rs_gmat(self, nu_s, temperature, vs=3, js=70, branches=(0,),
                    del_Tv=0.):
        """Resonant susceptibility based on G-matrix.

        This is implemented following :cite:`Koszykowski:85`.

        Parameters
        ----------
        nu_s : 1d array of floats
            Stokes frequencies (i.e., the calculation spectral domain).
        temperature : float
            Temperature in [K].
        vs : int, optional
            Total number of vibrational levels to be considered, by default 3.
        js : int, optional
            Total number of rotational levels to be considered, by default 70.
        branches : list of int, optional
            Branches to be considered, by default (0,) (i.e., only Q-branch).
        del_Tv : float, optional
            Absolute difference between vibrational and rotational temperature,
            by default 0.

        Returns
        -------
        1d array (complex)
            Theoretical resonant contributions (complex) based on G-matrix.
        """
        # Construct the v-branch-independent relaxation rate matrix
        gamma_mat = self.relax_mat(temperature, js)

        # For different v-branch combinations
        _js = np.arange(js)
        chi_rs = np.zeros_like(nu_s, dtype='complex128')
        for _branch in branches:
            for _v in np.arange(vs):
                # Calculate line positions
                nu_raman = self.ls_factors.line_pos(_v, _js, branch=_branch)

                # Construct the K_mat
                K_mat = np.diag(nu_raman) + gamma_mat*1j

                # Solve eigenvalue problem of K_mat
                eigvals, eigvec = np.linalg.eig(K_mat)
                eigvec_inv = np.linalg.inv(eigvec)

                # Compute the resonant intensity
                del_pop = self.ls_factors.pop_factor(temperature, _v, _js,
                                                     branch=_branch,
                                                     del_Tv=del_Tv)
                d = (self.trans_amp(_v, _js, branch=_branch))**0.5
                _term_l = d @ eigvec
                _term_r = eigvec_inv @ np.diag(del_pop) @ d
                _term = _term_l*_term_r

                for _j in _js:
                    _term_b = ((-nu_s + np.real(eigvals[_j]))**2
                               + np.imag(eigvals[_j])**2)
                    # A 1/2 factor is necessary to match the magnitude from
                    # isolated line assumption
                    chi_rs += 1/2*_term[_j]*np.conj(
                        -nu_s + eigvals[_j])/_term_b

        # A factor of c [cm/s] needs to be considered to convert cm^-1 to s^-1
        # by 2*pi*c
        return chi_rs/2/np.pi/self.C

    def chi_rs_isolated(self, nu_s, temperature, vs=3, js=70,
                        branches=(0, 2, -2), del_Tv=0):
        """Resonant susceptibility based on isolated line approximation.

        Parameters
        ----------
        nu_s : 1d array of floats
            Stokes frequencies (i.e., the calculation spectral domain).
        temperature : float
            Temperature in the probe volume.
        vs : int, optional
            Total number of vibrational levels to be considered, by default 3.
        js : int, optional
            Total number of rotational levels to be considered, by default 70.
        branches : list of int, optional
            Branches to be considered, by default (2, -2, 0).
        del_Tv : float, optional
            Absolute difference between vibrational and rotational temperature,
            by default 0.

        Returns
        -------
        1d array, 1d array
            Theoretical resonant contributions (complex) based on isolated line
            assumption.
        """
        # Construct the v-branch-independent relaxation rate matrix
        gamma_mat = self.relax_mat(temperature, js)
        Gamma_js = 2*abs(np.diag(gamma_mat))

        # Sum over all the transitions considered
        _js = np.arange(js)
        chi_rs = np.zeros_like(nu_s, dtype='complex128')
        for _branch in branches:
            for _v in np.arange(vs):
                # Calculate line positions, population factors, Gamma_j for all
                # js:
                nu_raman = self.ls_factors.line_pos(_v, _js, branch=_branch)
                del_pop = self.ls_factors.pop_factor(temperature, _v, _js,
                                                     branch=_branch,
                                                     del_Tv=del_Tv)
                d_squared = self.trans_amp(_v, _js, branch=_branch)

                # Calculate resonant contribution from each transition
                for _j in _js:
                    del_nu = nu_raman[_j] - nu_s
                    chi_rs += del_pop[_j]*d_squared[_j]*(
                        2*del_nu + Gamma_js[_j]*1j)/(
                            4*del_nu**2 + Gamma_js[_j]**2)

        # A factor of c [cm/s] needs to be considered to convert cm^-1 to s^-1
        # by 2*pi*c
        return chi_rs/2/np.pi/self.C

    def signal_as(self, temperature, x_mol=None, nu_s=None, pump_lw=None,
                  synth_mode=None, eq_func=eq_comp, **kwargs):
        r"""Anti-Stokes signal convoluted with a chosen laser lineshape.

        .. note::

            The Kataoka (or T-K) convolution :cite:`Kataoka:82, Teets:84` is a
            simplified form derived by :cite:`Farrow:85`. It is essentially
            an average of convolve(chi, pump)**2 and convolve(chi**2, pump),
            which takes care of the cross coherence effect. Excellent
            agreements are found when compared with Figs.13-15 in the NRC
            report :cite:`Parameswaran:88`.

        Parameters
        ----------
        temperature : float
            Temperature in the probe volume.
        x_mol : float, optional
            Mole fraction of probed molecule, by default None.
        nu_s : 1d array of float, optional
            Stokes frequencies (i.e., the calculation spectral domain), by
            default None.
        pump_lw : float, optional
            Pump laser linewdith (FWHM) in [:math:`\mathrm{cm}^{-1}`],
            by default None (i.e., no laser convolution is performed).
        synth_mode : dict, optional
            A dictionary containing the control parameters for creating the
            CARS spectrum, by default:

            pump_ls : 'Gaussian'
                Choose between pump laser lineshape between 'Gaussian' and
                'Lorentzian'.
            chi_rs : 'G-matrix'
                The method to compute resonant susceptibility: 'isolated' or
                'G-matrix' (collisional narrowing).
            convol : 'Kataoka' or 'K'
                Two ways to convolve the laser line with the resonant CARS
                susceptibilities:

                * 'Yuratich' or 'Y': One convolution, no cross coherence
                  effects accounted for :cite:`Yuratich:79`.
                * 'Kataoka' or 'K': Cross-coherence convolution, following
                  :cite:`Kataoka:82, Teets:84`.
                  This is necessary if pump linewidth is comparable to the
                  Raman linewidth (as is often the case at low T) and if
                  nonresonant contribution competes with resonant signal.
            doppler_effect : True
                Whether or not to take into account additional Doppler
                broadening.
            chem_eq : False
                Whether or not to assume chemical equilibrium. If True, an
                `eq_func` needs to be provided.
        eq_func : func, optional
                A function used to calculate local gas composition based on
                temperature and initial gas composition. Default is a simple
                template employing `cantera`.

        Other Parameters
        ----------------
        vs : int, optional
            Total number of vibrational levels to be considered, by default 3.
        js : int, optional
            Total number of rotational levels to be considered, by default 70.
        branches : list of int, optional
            Branches to be considered, by default (2, -2, 0).
        del_Tv : float, optional
            Absolute difference between vibrational and rotational temperature,
            by default 0.

        Returns
        -------
        1d array of floats, 1d array of floats
            Spectral locations in wavenumbers and theoretical CARS spectrum. If
            pump_lw is not specified, I_as needs to be multiplied by 1e-30 to
            retain its actual physical magnitude in
            [:math:`\mathrm{cm}^3/\mathrm{erg}`].
        """
        if synth_mode is None:
            self.synth_mode = {'pump_ls': 'Gaussian',
                               'chi_rs': 'G-matrix',
                               'convol': 'Kataoka',
                               'doppler_effect': False,
                               'chem_eq': False,
                               }
        else:
            self.synth_mode = synth_mode

        if x_mol is None:
            x_mol = self.init_comp[self.ls_factors.species]

        if self.synth_mode['chem_eq']:
            gas_comp = eq_func(temperature, self.pressure, self.init_comp)
            gas_comp = comp_normalize(gas_comp)
            N_species = gas_comp[
                self.ls_factors.species]*self.num_dens(temperature)
            chi_nrs = self.chi_nrs_est(temperature, local_comp=gas_comp)*1e-18
            self.local_comp = gas_comp
        else:
            # Species number density
            N_species = x_mol*self.num_dens(temperature)
            # Estimated nonresonant contribution based on given gas composition
            chi_nrs = self.chi_nrs_est(temperature)*1e-18

        # Specify the spectral domain
        if nu_s is None:
            # Ideally the grid should be much finer to correctly cover the
            # high-J transitions
            nu_s = np.linspace(2260, 2345, num=10000)

        # Compute the resonant contribution
        if self.synth_mode['chi_rs'] == 'isolated':
            chi_rs = self.chi_rs_isolated(nu_s, temperature, **kwargs)
        elif self.synth_mode['chi_rs'] == 'G-matrix':
            chi_rs = self.chi_rs_gmat(nu_s, temperature, **kwargs)
        else:
            raise ValueError("Unknown method. Only 'isolated' or 'G-matrix' "
                             "are available")

        # Compute the signal without pump laser convolution
        chi = (N_species*chi_rs + chi_nrs)*1e15

        # Consider Doppler broadening
        if self.synth_mode['doppler_effect']:
            _lw = self.ls_factors.doppler_lw(temperature)
            _ls = gaussian_line(nu_s, (nu_s[0]+nu_s[-1])/2, _lw)

            chi = np.convolve(chi, _ls, 'same')

        I_as = np.abs(chi)**2

        # Convolution with pump laser profile
        if pump_lw is not None:
            pump = []
            # Compute the pump laser line
            if self.synth_mode['pump_ls'] == 'Gaussian':
                pump = gaussian_line(nu_s, (nu_s[0]+nu_s[-1])/2, pump_lw)
            elif self.synth_mode['pump_ls'] == 'Lorentzian':
                pump = lorentz_line(nu_s, (nu_s[0]+nu_s[-1])/2, pump_lw)

            # Compute the convolution with Yuratich method (one convolution)
            if self.synth_mode['convol'] in ('Yuratich', 'Kataoka', 'Y', 'K'):
                I_as = np.convolve(I_as, pump, 'same')

            # Compute the convolution with Kataoka method (two convolutions)
            if self.synth_mode['convol'] in ('Kataoka', 'K'):
                chi_convol = np.convolve(chi, pump, 'same')
                # A correction factor is needed for the convolve(A)**2 v.s.
                # convolve(A**2)
                _corr = nu_s[1] - nu_s[0]
                I_as = 1/2*(I_as + _corr*np.abs(chi_convol)**2)

        return nu_s, I_as

ext: py

max_stars_repo_name: chuckedfromspace/carspy

max_issues_repo_name: chuckedfromspace/carspy

max_forks_repo_name: chuckedfromspace/carspy

reason: import numpy


text:
"""
Astra Viso demos.
"""

import astraviso as av
import numpy as np

def default():
    """
    Default Astra Viso demo. Creates a single image.

    Parameters
    ----------
    None

    Returns
    -------
    None
    """

    # Create StarCam and define motion
    cam = av.StarCam()
    cam.set_pointing_preset("kinematic", initial_quaternion=[0,0,0,1,],                            \
                                                                 initial_angular_rate=[0,0.025,0.1])

    # Setup star catalog
    cam.star_catalog.load_preset("random", 100000)

    # Create WorldObject and define motion
    cam.add_worldobject()
    cam.external_objects[0].set_position_preset("kinematic",                                       \
                 initial_position=np.array([0, 0, 1000]), initial_velocity=np.array([50, 25, 0]))
    cam.external_objects[0].set_vismag_preset("constant", vismag=2)

    # Create image and display
    av.imageutils.imshow(cam.snap(0, 1), [])

def sequence(count=2):
    """
    Astra Viso sequence demo.

    Parameters
    ----------
    count : int, optional
        Number of sequential images to generate. Default is 2.

    Returns
    -------
    None
    """

    # Create StarCam and define motion
    cam = av.StarCam()
    cam.set_pointing_preset("kinematic", initial_quaternion=[0,0,0,1,],                            \
                                                                  initial_angular_rate=[0.01, 0, 0])

    # Setup star catalog
    cam.star_catalog.load_preset("random", 100000)

    # Create WorldObject and define motion
    cam.add_worldobject()
    cam.external_objects[0].set_position_preset("kinematic",                                       \
                         initial_position=np.array([0, 0, 1]), initial_velocity=np.array([0, 0, 0]))
    cam.external_objects[0].set_vismag_preset("constant", vismag=-1)

    # Create image and display
    for idx in range(count):
        av.imageutils.imshow(cam.snap(idx, 1), [])

max_forks_repo_name: bradsease/test


text:
import networkx as nx

def getStat(G):
    print " \t\t\t STATISTICS ABOUT THE GRAPH:\n"
    print " - Number of nodes: \t", len(G.nodes())
    print " - Number of edges: \t", len(G.edges())
    print " - Number of isola: \t", len(nx.isolates(G))
    #print " - Number of trian: \t", nx.triangles(G)
    #print " - Clusturing coef: \t", nx.average_clustering(G)
    #print " - Transitive rato: \t", nx.transitivity(G)

reason: import networkx


text:
from wobbles.disc import Disc
import numpy as np

def compute_df_time_dependent(potential_extension_local, satellite_integration_time_list,
                              satellite_orbit_list, satellite_potential_list, velocity_dispersion_local,
                              rho_midplane=None, component_amplitude=None, verbose=False):

    """
    Does this seem right?
    :param potential_extension_local:
    :param orbit_integration_time_list:
    :param conversion_to_internal_time:
    :param satellite_orbit_list:
    :param satellite_potential_list:
    :param velocity_dispersion_local:
    :param rho_midplane:
    :return:
    """

    df_list, dj_list, force_list = [], [], []
    disc_instance = Disc(potential_extension_local)

    for satellite_integration_time in satellite_integration_time_list:

        t_end_sat = satellite_integration_time[-1]-satellite_integration_time[0]
        t_eval_orbits = np.linspace(0., t_end_sat, len(satellite_integration_time))

        df,  dj, f = compute_df(disc_instance, satellite_integration_time, satellite_orbit_list,
                                satellite_potential_list, velocity_dispersion_local, component_amplitude,
                                rho_midplane, t_eval_orbits, verbose)
        df_list.append(df)
        dj_list.append(dj)
        force_list.append(f)

    return df_list, dj_list, force_list

def compute_df(disc, t_eval_satellite, satellite_orbit_list, satellite_potential_list, velocity_dispersion_local,
                component_amplitude=None, component_densities=None, rho_midplane=None,
               t_eval_orbits=None, verbose=False, kwargs_distribution_function={}):

    """
    This function executes a certain workflow sequence: From a the orbit of a passing satellite, compute the
    perturbation to the action in the local phase space, and use it to derive a distribution function for the density/velocity
    in the z direction in the solar neighborhood.

    :param disc: An instance of Disc initialized with a local potential (see documentation in wobbles.wobbles)
    Note that this can be different from the galactic potential used to integrate the orbit of the passing satellite
    :param t_eval_satellite: The time over which to compute the perturbation from the satellite in internal units
    :param satellite_orbit_list: A list of perturbing satellite orbits
    :param satellite_potential_list: A list of satellite potentials corresponding to each orbit
    :param velocity_dispersion_local: The local velocity dispersion [km/sec]; if specified as a list, it corresponds to
    each component of a distribution function and must be the same length as component_amplitude
    :param component_amplitude: a list of amplitudes for each component in the disk; if specified must sum to one and be
    the same length as velocity dispersion local
    :param component_densities [M_solar/pc^3]: a list of densities corresponding to each velocity dispersion specified in velocity_dispersion_local
    This argument is ignored if
    A) velocity dispersion local is not a list, in which case the midplane density is computed as rho_midplane
    B) component amplitude is specified, in which case the component densities are computed as
    [amp_1 * rho_midplane, amp_2 * rho_midplane... ] where amp_1, amp_2... are specified in component_amplitude

    :param rho_midplane: the midplane density of the disk, needs to be specified for Isothermal potentials. For others it can
    be directly computed from the local potential (see rho_midplane method in potential_extension class)
    :param t_eval_orbits: the times when to evaluate the orbits of test particles in phase space; if None, reverts to the
    same times as those used to sample the satellite orbit
    :param verbose: makes print statements appear
    :param kwargs_distribution_function: additional keyword arguments for the DistributionFunction class

    :return: The distribution function, the perturbation to the action, and the force from the passing satellite
    """

    if t_eval_orbits is None:
        t_eval_orbits = t_eval_satellite

    if verbose:
        print('computing the force from ' + str(len(satellite_orbit_list)) + ' satellite orbits...')

    disc_phase_space_orbits = disc.orbits_in_phase_space(t_eval_orbits)

    force = disc.satellite_forces(t_eval_satellite, t_eval_orbits, satellite_orbit_list, satellite_potential_list,
                                  disc_phase_space_orbits, verbose)

    delta_J = disc.action_impulse(force, t_eval_orbits, satellite_orbit_list, satellite_potential_list,
                                 disc_phase_space_orbits)

    dF = disc.distribution_function(delta_J, velocity_dispersion_local, rho_midplane, component_amplitude, component_densities,
                                    verbose, kwargs_distribution_function)

    return dF, delta_J, force

max_stars_repo_name: dangilman/wobbles

max_issues_repo_name: dangilman/wobbles

max_forks_repo_name: dangilman/wobbles


text:
import ast
import os
import pickle
import sys

from PIL import Image as im
import numpy as np
from theanet.neuralnet import NeuralNet

from banti.iast_unicodes import LabelToUnicodeConverter
from banti.scaler import ScalerFactory
from banti.proglyph import ProGlyph, Space
from banti.processedpage import ProcessedPage
from scripts.tile import tile_raster_images

############################################# Arguments

if len(sys.argv) < 5:
    print("Usage:"
    "\n{0} neuralnet_params.pkl inpt.box/tiff scaler_params.scl codes.lbl "
    "\n\te.g:- {0} 0default.pkl sample_images/praasa.box "
    "scalings/relative48.scl labellings/alphacodes.lbl"
    "".format(sys.argv[0]))
    sys.exit()

nnet_prms_file_name = sys.argv[1]
input_file_name = sys.argv[2]
scaler_prms_file = sys.argv[3]
labelings_file_name = sys.argv[4]

############################################# Load Params

with open(scaler_prms_file, 'r') as sfp:
    scaler_prms = ast.literal_eval(sfp.read())

with open(nnet_prms_file_name, 'rb') as nnet_prms_file:
    nnet_prms = pickle.load(nnet_prms_file)

with open(labelings_file_name, encoding='utf-8') as labels_fp:
    labellings = ast.literal_eval(labels_fp.read())

# print(labellings)
chars = LabelToUnicodeConverter(labellings).onecode

############################################# Init Network
ProGlyph.scaler = ScalerFactory(scaler_prms)
procd_page = ProcessedPage(input_file_name)

nnet_prms['training_params']['BATCH_SZ'] = 1
ntwk = NeuralNet(**nnet_prms)
tester = ntwk.get_data_test_model(go_nuts=True)

############################################# Image saver
dir_name = os.path.basename(nnet_prms_file_name)[:-7] + '/'
if not os.path.exists(dir_name):
    os.makedirs(dir_name)
namer = (dir_name + '{:03d}_{}_{:02d}.png').format
print("Look for me in :", dir_name)

def saver(outs, ch, debug=True):
    saver.index += 1
    for i, out in enumerate(outs):
        global_normalize = False
        if out.ndim == 2:
            n_nodes = out.shape[1]
            w = n_nodes // int(np.sqrt(n_nodes))
            h = np.ceil(float(n_nodes) / w)
            extra = np.full((1, w*h-n_nodes), 0)
            out = np.concatenate((out, extra), 1).reshape((1, h, w))
        elif out.ndim == 4:
            out = out[0]
            if out.shape[-1] * out.shape[-2] < 65:
                global_normalize = True

        if debug:
            print("{:6.2f} {:6.2f} {:6.2f} {} GNM:{}".format(
                out.max(), out.mean(), out.min(), out.shape,
                global_normalize))

        im.fromarray(tile_raster_images(out, zm=2,
                                        make_white=True,
                                        global_normalize=global_normalize)
        ).save(namer(saver.index, chars[ch], i), compress_level=1)

    if debug:
        print()

saver.index = 0

############################################# Read glyphs & classify
print("Classifying...")
for line_pglyphs in procd_page.file_glyphs:
    for pglyph in line_pglyphs:
        if pglyph is Space:
            continue

        scaled_glp = pglyph.scaled
        img = scaled_glp.pix.astype('float32').reshape((1,)+scaled_glp.pix.shape)

        if ntwk.takes_aux():
            dtopbot = scaled_glp.dtop, scaled_glp.dbot
            aux_data = np.array([[dtopbot, dtopbot]], dtype='float32')
            logprobs_or_feats, preds, *layer_outs = tester(img, aux_data)
        else:
            logprobs_or_feats, preds, *layer_outs = tester(img)

        saver(layer_outs, np.argmax(logprobs_or_feats))
    print("Saved images of {} glyphs".format(saver.index))

print("Look for me in :", dir_name)

max_stars_repo_name: aimsravi/TELUGU-OCR

max_issues_repo_name: aimsravi/TELUGU-OCR

max_forks_repo_name: aimsravi/TELUGU-OCR


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

# # Desafio 4
#
# Neste desafio, vamos praticar um pouco sobre testes de hipóteses. Utilizaremos o _data set_ [2016 Olympics in Rio de Janeiro](https://www.kaggle.com/rio2016/olympic-games/), que contém dados sobre os atletas das Olimpíadas de 2016 no Rio de Janeiro.
#
# Esse _data set_ conta com informações gerais sobre 11538 atletas como nome, nacionalidade, altura, peso e esporte praticado. Estaremos especialmente interessados nas variáveis numéricas altura (`height`) e peso (`weight`). As análises feitas aqui são parte de uma Análise Exploratória de Dados (EDA).
#
# > Obs.: Por favor, não modifique o nome das funções de resposta.

# ## _Setup_ geral

# In[1]:

import pandas as pd
import matplotlib.pyplot as plt
import numpy as np
import scipy.stats as sct
import seaborn as sns
import statsmodels.api as sm

# In[2]:

#%matplotlib inline
from IPython.core.pylabtools import figsize

figsize(12, 8)

sns.set()

# In[3]:

athletes = pd.read_csv("athletes.csv")

# In[4]:

def get_sample(df, col_name, n=100, seed=42):
    """Get a sample from a column of a dataframe.

    It drops any numpy.nan entries before sampling. The sampling
    is performed without replacement.

    Example of numpydoc for those who haven't seen yet.

    Parameters
    ----------
    df : pandas.DataFrame
        Source dataframe.
    col_name : str
        Name of the column to be sampled.
    n : int
        Sample size. Default is 100.
    seed : int
        Random seed. Default is 42.

    Returns
    -------
    pandas.Series
        Sample of size n from dataframe's column.
    """
    np.random.seed(seed)

    random_idx = np.random.choice(df[col_name].dropna().index, size=n, replace=False)

    return df.loc[random_idx, col_name]

# ## Inicia sua análise a partir daqui

# In[5]:

h_samp = get_sample(athletes, 'height', 3000)
w_samp = get_sample(athletes, 'weight', 3000)
alpha = 0.05
w_log = np.log(w_samp)

# ## Questão 1
#
# Considerando uma amostra de tamanho 3000 da coluna `height` obtida com a função `get_sample()`, execute o teste de normalidade de Shapiro-Wilk com a função `scipy.stats.shapiro()`. Podemos afirmar que as alturas são normalmente distribuídas com base nesse teste (ao nível de significância de 5%)? Responda com um boolean (`True` ou `False`).

# In[6]:

def q1():
    shapiro = sct.shapiro(h_samp)
    print(shapiro[1])
    return bool(shapiro[1] > alpha)

# In[7]:

q1()

# __Para refletir__:
#
# * Plote o histograma dessa variável (com, por exemplo, `bins=25`). A forma do gráfico e o resultado do teste são condizentes? Por que?
# * Plote o qq-plot para essa variável e a analise.
# * Existe algum nível de significância razoável que nos dê outro resultado no teste? (Não faça isso na prática. Isso é chamado _p-value hacking_, e não é legal).

# In[8]:

h_samp.plot(kind='hist', bins=25);
# Resposta: Não, pois a forma do gráfico me parece normal, mas o resultado do teste, falso, indica
# a evidência de que os dados não vêm de uma distribuição normal.

# In[9]:

sm.qqplot(h_samp, fit=True, line='45');
# Resposta: Não há, pois conforme aula, os valores típicos de α/alpha são 0.025, 0.05 e 0.10 e o P neste caso é bem menor.

# ## Questão 2
#
# Repita o mesmo procedimento acima, mas agora utilizando o teste de normalidade de Jarque-Bera através da função `scipy.stats.jarque_bera()`. Agora podemos afirmar que as alturas são normalmente distribuídas (ao nível de significância de 5%)? Responda com um boolean (`True` ou `False`).

# In[10]:

def q2():
    jBera = sct.jarque_bera(h_samp)
    print(jBera[1])
    return bool(jBera[1] > alpha)

# In[11]:

q2()

# __Para refletir__:
#
# * Esse resultado faz sentido?

# In[12]:

# Resposta: creio que sim, pois novamente indica que a amostra não é normal, pois se fosse, P deveria estar mais próximo de 0.5.

# ## Questão 3
#
# Considerando agora uma amostra de tamanho 3000 da coluna `weight` obtida com a função `get_sample()`. Faça o teste de normalidade de D'Agostino-Pearson utilizando a função `scipy.stats.normaltest()`. Podemos afirmar que os pesos vêm de uma distribuição normal ao nível de significância de 5%? Responda com um boolean (`True` ou `False`).

# In[13]:

def q3():
    pearson = sct.normaltest(w_samp)
    print(pearson[1])
    return bool(pearson[1] > alpha)

# In[14]:

q3()

# __Para refletir__:
#
# * Plote o histograma dessa variável (com, por exemplo, `bins=25`). A forma do gráfico e o resultado do teste são condizentes? Por que?
# * Um _box plot_ também poderia ajudar a entender a resposta.

# In[15]:

w_samp.plot(kind='hist', bins=25);
# Resposta: agora sim, pois o gráfico não indica uma curva normal e o teste fez o mesmo.

# ## Questão 4
#
# Realize uma transformação logarítmica na amostra de `weight` da questão 3 e repita o mesmo procedimento. Podemos afirmar a normalidade da variável transformada ao nível de significância de 5%? Responda com um boolean (`True` ou `False`).

# In[16]:

def q4():
    pearson = sct.normaltest(w_log)
    print(pearson[1])
    return bool(pearson[1] > alpha)

# In[17]:

q4()

# __Para refletir__:
#
# * Plote o histograma dessa variável (com, por exemplo, `bins=25`). A forma do gráfico e o resultado do teste são condizentes? Por que?
# * Você esperava um resultado diferente agora?

# In[18]:

w_log.plot(kind='hist', bins=25);
# Resposta: não, pois novamente a curva é normal.

# > __Para as questão 5 6 e 7 a seguir considere todos testes efetuados ao nível de significância de 5%__.

# In[19]:

athletes.head()

# In[20]:

h_bra = athletes.query('nationality == "BRA"')['height'].dropna()
h_usa = athletes.query('nationality == "USA"')['height'].dropna()
h_can = athletes.query('nationality == "CAN"')['height'].dropna()
print(h_bra.size, h_usa.size, h_can.size)

# ## Questão 5
#
# Obtenha todos atletas brasileiros, norte-americanos e canadenses em `DataFrame`s chamados `bra`, `usa` e `can`,respectivamente. Realize um teste de hipóteses para comparação das médias das alturas (`height`) para amostras independentes e variâncias diferentes com a função `scipy.stats.ttest_ind()` entre `bra` e `usa`. Podemos afirmar que as médias são estatisticamente iguais? Responda com um boolean (`True` ou `False`).

# In[21]:

def q5():
    bra_usa = sct.ttest_ind(h_bra, h_usa, equal_var=False)
    print(bra_usa[1])
    return bool(bra_usa[1] > alpha)

# In[22]:

q5()

# ## Questão 6
#
# Repita o procedimento da questão 5, mas agora entre as alturas de `bra` e `can`. Podemos afimar agora que as médias são estatisticamente iguais? Reponda com um boolean (`True` ou `False`).

# In[23]:

def q6():
    bra_can = sct.ttest_ind(h_bra, h_can, equal_var=False)
    print(bra_can[1])
    return bool(bra_can[1] > alpha)

# In[24]:

q6()

# ## Questão 7
#
# Repita o procedimento da questão 6, mas agora entre as alturas de `usa` e `can`. Qual o valor do p-valor retornado? Responda como um único escalar arredondado para oito casas decimais.

# In[25]:

def q7():
    usa_can = sct.ttest_ind(h_usa, h_can, equal_var=False)
    return float(round(usa_can[1], 8))

# In[26]:

q7()

# __Para refletir__:
#
# * O resultado faz sentido?
# * Você consegue interpretar esse p-valor?
# * Você consegue chegar a esse valor de p-valor a partir da variável de estatística?

reason: import numpy,import scipy,import statsmodels


text:
try:
	from .model_interface import model_zoo
except:
	from model_interface import model_zoo

import tensorflow as tf
import numpy as np
from utils.bert import bert_utils
from loss import loss_utils, triplet_loss_utils

from model_io import model_io
from task_module import classifier
import tensorflow as tf
from metric import tf_metrics

from optimizer import distributed_optimizer as optimizer
from model_io import model_io
from utils.simclr import simclr_utils

def get_labels_of_similarity(query_input_ids, anchor_query_ids):
	idxs_1 = tf.expand_dims(query_input_ids, axis=1) # batch 1 seq
	idxs_2 = tf.expand_dims(anchor_query_ids, axis=0) # 1 batch seq
	# batch x batch x seq
	labels = tf.cast(tf.not_equal(idxs_1, idxs_2), tf.float32) # not equal:1, equal:0
	equal_num = tf.reduce_sum(labels, axis=-1) # [batch, batch]
	not_equal_label = tf.cast(tf.not_equal(equal_num, 0), tf.float32)
	not_equal_label_shape = bert_utils.get_shape_list(not_equal_label, expected_rank=[2,3])
	not_equal_label *= tf.cast(1 - tf.eye(not_equal_label_shape[0]), tf.float32)
	equal_label = (1 - not_equal_label) - tf.eye(not_equal_label_shape[0])
	return equal_label, not_equal_label

def get_finised_pos_v1(token_seq, finished_index, max_length):
	token_seq = tf.cast(token_seq, tf.int32)
	seq_shape = bert_utils.get_shape_list(token_seq, expected_rank=[2,3])
	match_indices = tf.where(                          # [[5, 5, 2, 5, 4],
	tf.equal(finished_index, token_seq),                              #  [0, 5, 2, 3, 5],
		x=tf.range(seq_shape[1]) * tf.ones_like(token_seq),  #  [5, 1, 5, 5, 5]]
		y=(seq_shape[1])*tf.ones_like(token_seq))

	finished_pos = tf.reduce_min(match_indices, axis=1)
	mask = tf.cast(tf.one_hot(finished_pos, max_length), tf.float32) # [batch, max_length]

	modified_token_seq = tf.cast(token_seq, tf.float32) - float(finished_index) * mask + mask * 1.0

	return tf.cast(modified_token_seq, tf.int32)

def model_fn_builder(
					model_config,
					num_labels,
					init_checkpoint,
					model_reuse=None,
					load_pretrained=True,
					model_io_config={},
					opt_config={},
					exclude_scope="",
					not_storage_params=[],
					target="a",
					label_lst=None,
					output_type="sess",
					**kargs):

	def model_fn(features, labels, mode):

		model_api = model_zoo(model_config)

		# if kargs.get('trf_input', False):

		# 	input_shape_list = bert_utils.get_shape_list(features['input_ids'], expected_rank=2)
		# 	batch_size = input_shape_list[0]
		# 	seq_length = input_shape_list[1]

		# 	input_ids = get_finised_pos_v1(features['input_ids'], 102, seq_length)
		# 	features['input_ids'] = input_ids

		# 	tf.logging.info("**** trf input modification for qa and sentence pair **** ")

		# features['segment_ids'] = 0 * features['segment_ids']
		model = model_api(model_config, features, labels,
							mode, target, reuse=model_reuse, **kargs)

		label_ids = features["label_ids"]

		if mode == tf.estimator.ModeKeys.TRAIN:
			dropout_prob = model_config.get("dropout_prob", 0.1)
			is_training = True
		else:
			dropout_prob = 0.0
			is_training = False

		if model_io_config.fix_lm == True:
			scope = model_config.scope + "_finetuning"
		else:
			scope = model_config.scope

		with tf.variable_scope(scope, reuse=model_reuse):
			(loss,
				per_example_loss,
				logits) = classifier.classifier(model_config,
											model.get_pooled_output(),
											num_labels,
											label_ids,
											dropout_prob)

		if not kargs.get('use_tpu'):
			tf.summary.scalar("classifier_loss", loss)

		cpc_flag = model_config.get('apply_cpc', 'none')
		if model_config.get("label_type", "single_label") == "multi_label":
			pos_true_mask, neg_true_mask = get_labels_of_similarity(label_ids, label_ids)
		else:
			label_shape = bert_utils.get_shape_list(label_ids, expected_rank=[1, 2,3])
			neg_true_mask = tf.cast(triplet_loss_utils._get_anchor_negative_triplet_mask(label_ids), tf.float32)
			pos_true_mask = (1.0 - neg_true_mask)*(1.0-tf.eye(label_shape[0], dtype=tf.float32))
		if cpc_flag == 'apply':
			loss_mask = tf.minimum(tf.reduce_sum(pos_true_mask, axis=-1), 1)
			feat = model.get_pooled_output()

			if kargs.get("apply_feat_projection", False):
				with tf.variable_scope(scope+"/head_proj", reuse=tf.AUTO_REUSE):
					feat = simclr_utils.projection_head(feat,
											is_training,
											head_proj_dim=128,
											num_nlh_layers=1,
											head_proj_mode='nonlinear',
											name='head_contrastive')
				feat = tf.nn.l2_normalize(feat+1e-20, axis=-1)
				cosine_score = tf.matmul(feat, tf.transpose(feat)) / 0.1
				tf.logging.info("****** apply simclr projection and. l2 normalize *******")
			else:
				feat = tf.nn.l2_normalize(feat+1e-20, axis=-1)
				cosine_score = tf.matmul(feat, tf.transpose(feat)) / 0.1
				tf.logging.info("****** apply raw feat *******")
			cosine_score_neg = neg_true_mask * cosine_score
			cosine_score_pos = -pos_true_mask * cosine_score

			y_pred_neg = cosine_score_neg - (1 - neg_true_mask) * 1e12
			y_pred_pos = cosine_score_pos - (1 - pos_true_mask) * 1e12

			# add circle-loss without margin and scale-factor
			joint_neg_loss = tf.reduce_logsumexp(y_pred_neg, axis=-1)
			joint_pos_loss = tf.reduce_logsumexp(y_pred_pos, axis=-1)
			cpc_loss = tf.reduce_sum(tf.nn.softplus(joint_neg_loss+joint_pos_loss)*loss_mask)/tf.reduce_sum(loss_mask+1e-10)
			loss += cpc_loss
			if not kargs.get('use_tpu'):
				tf.summary.scalar("cpc_loss", cpc_loss)
				tf.summary.scalar("cpc_mask", tf.reduce_mean(tf.reduce_sum(loss_mask, axis=-1)))
			tf.logging.info("****** apply cpc-loss *******")

		if kargs.get('apply_gp', False):
			gp_loss = loss_utils.gradient_penalty_loss(loss, model.get_embedding_table(),
											epsilon=1.0)
			loss += gp_loss
			tf.logging.info("****** apply gradient penalty *******")

		model_io_fn = model_io.ModelIO(model_io_config)

		tvars = model_io_fn.get_params(model_config.scope,
										not_storage_params=not_storage_params)

		if model_config.scope in ['textcnn']:
			emb_tvars = model_io_fn.get_params("textcnn_token_embedding",
										not_storage_params=not_storage_params)
			tvars.extend(emb_tvars)
		try:
			params_size = model_io_fn.count_params(model_config.scope)
			print("==total params==", params_size)
		except:
			print("==not count params==")
		print(tvars)
		if load_pretrained == "yes":
			model_io_fn.load_pretrained(tvars,
										init_checkpoint,
										exclude_scope=exclude_scope)

		if mode == tf.estimator.ModeKeys.TRAIN:

			optimizer_fn = optimizer.Optimizer(opt_config)

			model_io_fn.print_params(tvars, string=", trainable params")
			update_ops = tf.get_collection(tf.GraphKeys.UPDATE_OPS)
			print("==update_ops==", update_ops)
			with tf.control_dependencies(update_ops):
				train_op = optimizer_fn.get_train_op(loss, tvars,
								opt_config.init_lr,
								opt_config.num_train_steps,
								**kargs)

				model_io_fn.set_saver()

				if kargs.get("task_index", 1) == 0 and kargs.get("run_config", None):
					training_hooks = []
				elif kargs.get("task_index", 1) == 0:
					model_io_fn.get_hooks(kargs.get("checkpoint_dir", None),
														kargs.get("num_storage_steps", 1000))

					training_hooks = model_io_fn.checkpoint_hook
				else:
					training_hooks = []

				if len(optimizer_fn.distributed_hooks) >= 1:
					training_hooks.extend(optimizer_fn.distributed_hooks)
				print(training_hooks, "==training_hooks==", "==task_index==", kargs.get("task_index", 1))

				estimator_spec = tf.estimator.EstimatorSpec(mode=mode,
								loss=loss, train_op=train_op,
								training_hooks=training_hooks)
				print(tf.global_variables(), "==global_variables==")
				if output_type == "sess":
					return {
						"train":{
										"loss":loss,
										"logits":logits,
										"train_op":train_op
									},
						"hooks":training_hooks
					}
				elif output_type == "estimator":
					return estimator_spec

		elif mode == tf.estimator.ModeKeys.PREDICT:
			# if model_config.get('label_type', 'single_label') == 'single_label':
			# 	print(logits.get_shape(), "===logits shape===")
			# 	pred_label = tf.argmax(logits, axis=-1, output_type=tf.int32)
			# 	prob = tf.nn.softmax(logits)
			# 	max_prob = tf.reduce_max(prob, axis=-1)

			# 	estimator_spec = tf.estimator.EstimatorSpec(
			# 							mode=mode,
			# 							predictions={
			# 										'pred_label':pred_label,
			# 										"max_prob":max_prob
			# 							},
			# 							export_outputs={
			# 								"output":tf.estimator.export.PredictOutput(
			# 											{
			# 												'pred_label':pred_label,
			# 												"max_prob":max_prob
			# 											}
			# 										)
			# 							}
			# 				)
			if model_config.get('label_type', 'single_label') == 'multi_label':
				print("==apply multi_label==")
				if model_config.get("loss", "entropy") == "entropy":
					prob = tf.nn.sigmoid(logits)
				elif model_config.get("loss", "entropy") == "circle_loss":
					prob = tf.nn.sigmoid(logits)
				estimator_spec = tf.estimator.EstimatorSpec(
										mode=mode,
										predictions={
													'pred_label':prob,
													"max_prob":prob
										},
										export_outputs={
											"output":tf.estimator.export.PredictOutput(
														{
															'pred_label':prob,
															"max_prob":prob
														}
													)
										}
							)
			elif model_config.get('label_type', 'single_label') == "single_label":
				print("==apply multi_label==")
				prob = tf.nn.softmax(logits)
				estimator_spec = tf.estimator.EstimatorSpec(
										mode=mode,
										predictions={
													'pred_label':prob,
													"max_prob":prob
										},
										export_outputs={
											"output":tf.estimator.export.PredictOutput(
														{
															'pred_label':prob,
															"max_prob":prob
														}
													)
										}
							)
			return estimator_spec

		elif mode == tf.estimator.ModeKeys.EVAL:
			def metric_fn(per_example_loss,
						logits,
						label_ids):
				"""Computes the loss and accuracy of the model."""
				sentence_log_probs = tf.reshape(
					logits, [-1, logits.shape[-1]])
				sentence_predictions = tf.argmax(
					logits, axis=-1, output_type=tf.int32)
				sentence_labels = tf.reshape(label_ids, [-1])
				sentence_accuracy = tf.metrics.accuracy(
					labels=label_ids, predictions=sentence_predictions)
				sentence_mean_loss = tf.metrics.mean(
					values=per_example_loss)
				sentence_f = tf_metrics.f1(label_ids,
										sentence_predictions,
										num_labels,
										label_lst, average="macro")

				eval_metric_ops = {
									"f1": sentence_f,
									"acc":sentence_accuracy
								}

				return eval_metric_ops

			if output_type == "sess":
				return {
					"eval":{
							"per_example_loss":per_example_loss,
							"logits":logits,
							"loss":tf.reduce_mean(per_example_loss),
							"feature":model.get_pooled_output()
						}
				}
			elif output_type == "estimator":
				eval_metric_ops = metric_fn(
							per_example_loss,
							logits,
							label_ids)

				estimator_spec = tf.estimator.EstimatorSpec(mode=mode,
								loss=loss,
								eval_metric_ops=eval_metric_ops)
				return estimator_spec
		else:
			raise NotImplementedError()
	return model_fn

max_stars_repo_name: yyht/bert

max_issues_repo_name: yyht/bert

max_forks_repo_name: yyht/bert


text:
#original author -- Koki Shirota

import cv2
import numpy as np

cap = cv2.VideoCapture(0)
cap.set(3, 1280)
cap.set(4, 720)

while(1):
    ret, frame = cap.read()
    hsv = cv2.cvtColor(frame, cv2.COLOR_BGR2HSV)

    lower_blue = np.array([100, 50, 50])
    upper_blue = np.array([150, 255, 255])
    img_mask_blue = cv2.inRange(hsv, lower_blue, upper_blue)
    img_color_blue = cv2.bitwise_and(frame, frame, mask=img_mask_blue)

    cv2.imshow("SHOW COLOR DEFAULT", frame)
    cv2.imshow("SHOW COLOR BLUE", img_color_blue)

    k = cv2.waitKey(1)
    if k == 27:
        break

cap.release()
cv2.destroyAllWindows()


text:
#!/usr/bin/env python

"""
Perform data manipulation tasks in project workflow
"""

import os
try:
    import urllib.request as urllib2
except ImportError:
    import urllib2

import boto
import numpy as np
import pandas as pd
import boto.s3
from boto.s3.key import Key
from scipy.stats import randint as sp_randint
from sklearn.model_selection import RandomizedSearchCV
from sklearn.ensemble import RandomForestRegressor

__author__ = "Pearl Philip"
__credits__ = "David Beck"
__license__ = "BSD 3-Clause License"
__maintainer__ = "Pearl Philip"
__email__ = "pphilip@uw.edu"
__status__ = "Development"

def create_notation_dataframe(filename):
    """
    Returning Pandas dataframe of sample ID and molecular notation.
    :param filename: File object containing molecular notation indexed by sample ID
    :return: Dataframe of molecular notation indexed by sample ID.
    """
    df = []
    for line in filename:
        # Splits the line into it's key and molecular string
        words = line.split()
        z = [int(words[0]), words[1]]
        df.append(z)
    df = pd.DataFrame(df)
    df.columns = ['CID', 'SMILES']
    df.sort_values(by='CID', inplace=True)
    return df

def create_activity_dataframe(dataframe):
    """
    Performing useful transformations on the acquired data for use in subsequent algorithm.
    :param dataframe: Dataframe downloaded from NCBI database.
    :return: df: Cleaned and sorted dataframe.
    """

    # Eliminates first five text rows of csv
    for j in range(5):
        df = dataframe.drop(j, axis=0)
    df = df.drop(['PUBCHEM_ACTIVITY_URL', 'PUBCHEM_RESULT_TAG',
                  'PUBCHEM_ACTIVITY_SCORE', 'PUBCHEM_SID',
                  'PUBCHEM_ASSAYDATA_COMMENT', 'Potency',
                  'Efficacy', 'Analysis Comment',
                  'Curve_Description', 'Fit_LogAC50',
                  'Fit_HillSlope', 'Fit_R2', 'Fit_ZeroActivity',
                  'Fit_CurveClass', 'Excluded_Points', 'Compound QC',
                  'Max_Response', 'Phenotype', 'Activity at 0.457 uM',
                  'Activity at 2.290 uM', 'Activity at 11.40 uM',
                  'Activity at 57.10 uM', 'PUBCHEM_ACTIVITY_OUTCOME',
                  'Fit_InfiniteActivity'], axis=1)
    df.rename(columns={'PUBCHEM_CID': 'CID'}, inplace=True)

    # Eliminates duplicate compound rows
    df['dupes'] = df.duplicated('CID')
    df = df[df['dupes'] == 0].drop(['dupes'], axis=1)
    df = df.sort_values(by='CID')
    return df

def upload_to_s3(aws_access_key_id, aws_secret_access_key, file_to_s3, bucket, key, callback=None, md5=None,
                 reduced_redundancy=False, content_type=None):
    """
    Uploads the given file to the AWS S3 bucket and key specified.
    :param aws_access_key_id: First part of AWS access key.
    :param aws_secret_access_key: Second part of AWS access key.
    :param file_to_s3: File object to be uploaded.
    :param bucket: S3 bucket name as string.
    :param key: Name attribute of the file object to be uploaded.
    :param callback: Function accepts two integer parameters, the first representing the number of bytes that have been
    successfully transmitted to S3 and the second representing the size of the to be transmitted object. Returns
    boolean indicating success/failure of upload.
    :param md5: MD5 checksum value to verify the integrity of the object.
    :param reduced_redundancy: S3 option that enables customers to reduce their costs
    by storing noncritical, reproducible data at lower levels of redundancy than S3's standard storage.
    :param content_type: Set the type of content in file object.
    :return: Boolean indicating success of upload.
    """
    try:
        size = os.fstat(file_to_s3.fileno()).st_size
    except:
        # Not all file objects implement fileno(), so we fall back on this
        file_to_s3.seek(0, os.SEEK_END)
        size = file_to_s3.tell()

    conn = boto.connect_s3(aws_access_key_id, aws_secret_access_key)
    bucket = conn.get_bucket(bucket, validate=True)
    k = Key(bucket)
    k.key = key
    if content_type:
        k.set_metadata('Content-Type', content_type)
    sent = k.set_contents_from_file(file_to_s3, cb=callback, md5=md5,
                                    reduced_redundancy=reduced_redundancy, rewind=True)

    # Rewind for later use
    file_to_s3.seek(0)

    if sent == size:
        return True
    return False

def join_dataframes():
    """
    Joining the dataframes of existing descriptor files from their urls into a single dataframe.
    :return: Dataframe after join over key column.
    """

    url_list = ['https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_constitution.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_con.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_kappa.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_estate.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_basak.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_property.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_charge.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_moe.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_burden.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_geary.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_moran.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_topology.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_geometric.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_cpsa.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_rdf.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_morse.csv',
                'https://s3-us-west-2.amazonaws.com/pphilip-usp-inhibition/data/df_whim.csv'
                ]

    url_exist_list = []
    for url in url_list:
        try:
            r = urllib2.urlopen(url)
        except urllib2.URLError as e:
            r = e
        if r.code < 400:
            url_exist_list.append(url)

    i = 0
    df = [0] * len(url_exist_list)
    for url in url_exist_list:
        df[i] = pd.read_csv(url)

        df[i].drop(df[i].columns[0], axis=1, inplace=True)
        df[i].reset_index(drop=True, inplace=True)
        i += 1

    joined_df = df[0]
    for i in df[1:]:
        joined_df = joined_df.join(i)
    return joined_df

def choose_features(x_train, y_train, x_test, column_names):
    """
    Selecting the features of high importance to reduce feature space.
    :param x_train: Training set of features.
    :param x_test: Test set of features.
    :param y_train: Training target values
    :param column_names: Names of columns in x
    """
    # Random forest feature importance
    clf = RandomForestRegressor(n_jobs=-1, random_state=1, n_estimators=10)
    clf.fit(x_train, y_train.ravel())
    feature_importance = clf.feature_importances_
    scores_table = pd.DataFrame({'feature': column_names, 'scores':
                                 feature_importance}).sort_values(by=['scores'], ascending=False)
    scores = scores_table['scores'].tolist()
    n_features = [25, 50, 75, 100, 150, 200, 250, 300]
    for n in n_features:
        feature_scores = scores_table['feature'].tolist()
        selected_features = feature_scores[:n]
        x_train = pd.DataFrame(x_train, columns=column_names)
        desired_x_train = x_train[selected_features]
        x_test = pd.DataFrame(x_test, columns=column_names)
        desired_x_test = x_test[selected_features]

        desired_x_train.to_csv('../data/x_train_postprocessing_rfr_%d.csv' % n)
        desired_x_test.to_csv('../data/x_test_postprocessing_rfr_%d.csv' % n)
    pd.DataFrame(scores).to_csv('../data/feature_scores_rfr.csv')

    return

def change_nan_infinite(dataframe):
    """
    Replacing NaN and infinite values from the dataframe with zeros.
    :param dataframe: Dataframe containing NaN and infinite values.
    :return data: Data with no NaN or infinite values.
    """

    dataframe.replace([np.inf, -np.inf], np.nan, inplace=True)
    data = dataframe.fillna(0)

    return data

max_stars_repo_name: BeckResearchLab/USP-inhibition

max_issues_repo_name: BeckResearchLab/USP-inhibition

max_forks_repo_name: BeckResearchLab/USP-inhibition

reason: import numpy,from scipy