python类exp()的实例源码

def devi(yy, eta):
    deveta = yy*eta - scipy.exp(eta)
    devy = yy*scipy.log(yy) - yy
    devy[yy == 0] = 0
    result = 2*(devy - deveta)
    return(result)

def rbf_kernel_pca(X, gamma, n_components):
    sq_dists = pdist(X, 'sqeuclidean')
    mat_sq_dists = squareform(sq_dists)
    K = exp(-gamma * mat_sq_dists)

    N = K.shape[0]
    one_n = np.ones((N, N)) / N
    K = K - one_n.dot(K) - K.dot(one_n) + one_n.dot(K).dot(one_n)

    eigenvalues, eigenvectors = eigh(K)

    alphas = np.column_stack((
        eigenvectors[:, -i] for i in range(1, n_components+1)
    ))
    lambdas = [eigenvalues[-i] for i in range(1, n_components+1)]

    return alphas, lambdas

def KMV_f(E,D,T,r,sigmaE):
    n=10000
    m=2000
    diffOld=1e6     # a very big number
    for i in sp.arange(1,10):
        for j in sp.arange(1,m):
            A=E+D/2+i*D/n
            sigmaA=0.05+j*(1.0-0.001)/m
            d1 = (log(A/D)+(r+sigmaA*sigmaA/2.)*T)/(sigmaA*sqrt(T))
            d2 = d1-sigmaA*sqrt(T)
            diff4E= (A*N(d1)-D*exp(-r*T)*N(d2)-E)/A  # scale by assets
            diff4A= A/E*N(d1)*sigmaA-sigmaE          # a small number already
            diffNew=abs(diff4E)+abs(diff4A)
            if diffNew<diffOld:
               diffOld=diffNew
               output=(round(A,2),round(sigmaA,4),round(diffNew,5))
    return output
#

def implied_vol_call_min(S,X,T,r,c): 
    from scipy import log,exp,sqrt,stats 
    implied_vol=1.0
    min_value=1000
    for i in range(10000): 
        sigma=0.0001*(i+1)
        d1=(log(S/X)+(r+sigma*sigma/2.)*T)/(sigma*sqrt(T)) 
        d2 = d1-sigma*sqrt(T)
        c2=S*stats.norm.cdf(d1)-X*exp(-r*T)*stats.norm.cdf(d2) 
        abs_diff=abs(c2-c)
        if abs_diff<min_value: 
            min_value=abs_diff 
            implied_vol=sigma 
            k=i
    return implied_vol


# Step 3: get call option data

def implied_vol_put_min(S,X,T,r,p):
    implied_vol=1.0
    min_value=100.0
    for i in xrange(1,10000): 
        sigma=0.0001*(i+1)
        d1=(log(S/X)+(r+sigma*sigma/2.)*T)/(sigma*sqrt(T)) 
        d2 = d1-sigma*sqrt(T)
        put=X*exp(-r*T)*stats.norm.cdf(-d2)-S*stats.norm.cdf(-d1) 
        abs_diff=abs(put-p)
        if abs_diff<min_value: 
            min_value=abs_diff 
            implied_vol=sigma 
            k=i
        put_out=put
    print 'k, implied_vol, put, abs_diff' 
    return k,implied_vol, put_out,min_value

def up_and_out_call(s0,x,T,r,sigma,n_simulation,barrier):
    n_steps=100. 
    dt=T/n_steps 
    total=0 
    for j in sp.arange(0, n_simulation): 
        sT=s0 
        out=False
        for i in range(0,int(n_steps)): 
            e=sp.random.normal() 
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT>barrier: 
               out=True 
        if out==False: 
            total+=bsCall(s0,x,T,r,sigma) 
    return total/n_simulation 
#

def up_and_out_call(s0,x,T,r,sigma,n_simulation,barrier):
    n_steps=100. 
    dt=T/n_steps 
    total=0 
    for j in sp.arange(0, n_simulation): 
        sT=s0 
        out=False
        for i in range(0,int(n_steps)): 
            e=sp.random.normal() 
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT>barrier: 
               out=True 
        if out==False: 
            total+=bsCall(s0,x,T,r,sigma) 
    return total/n_simulation 
#

def upCall(s,x,T,r,sigma,nSimulation,barrier):
    import scipy as sp
    import p4f 
    n_steps=100
    dt=T/n_steps 
    inTotal=0 
    outTotal=0
    for j in range(0, nSimulation): 
        sT=s
        inStatus=False 
        outStatus=True
        for i in range(0,int(n_steps)):
            e=sp.random.normal()
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT>barrier:
                outStatus=False 
                inStatus=True
        if outStatus==True:
            outTotal+=p4f.bs_call(s,x,T,r,sigma) 
        else:
            inTotal+=p4f.bs_call(s,x,T,r,sigma) 
    return outTotal/nSimulation, inTotal/nSimulation
#

def lookback_min_price_as_strike(s,T,r,sigma,n_simulation): 
    n_steps=100
    dt=T/n_steps
    total=0
    for j in range(n_simulation): 
        min_price=100000.   # a very big number 
        sT=s
        for i in range(int(n_steps)): 
            e=sp.random.normal()
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT<min_price:
                min_price=sT
                #print 'j=',j,'i=',i,'total=',total 
                total+=p4f.bs_call(s,min_price,T,r,sigma)
    return total/n_simulation

#

def up_and_out_call(s0,x,T,r,sigma,n_simulation,barrier):
    n_steps=100. 
    dt=T/n_steps 
    total=0 
    for j in sp.arange(0, n_simulation): 
        sT=s0 
        out=False
        for i in range(0,int(n_steps)): 
            e=sp.random.normal() 
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT>barrier: 
               out=True 
        if out==False: 
            total+=bsCall(s0,x,T,r,sigma) 
    return total/n_simulation 
#

def up_and_out_call(s0,x,T,r,sigma,n_simulation,barrier):
    n_steps=100.
    dt=T/n_steps 
    total=0
    for j in range(0, n_simulation): 
        sT=s0
        out=False
        for i in range(0,int(n_steps)): 
            e=sp.random.normal()
            sT*=sp.exp((r-0.5*sigma*sigma)*dt+sigma*e*sp.sqrt(dt)) 
            if sT>barrier:
                out=True 
    if out==False:
       total+=p4f.bs_call(s0,x,T,r,sigma) 
    return total/n_simulation
#

def phaseFirstOrder(self, fromPos, toPos, phi1, degree = True, noChangeOnResonance=False, pivot = 0, scale= "Hz"):
        """This should only be applied to Fourier transformed spectral data
        It will lead to a zero phase shift at the upper end of the spectrum and to a
        phase shift of phi1 at the lower end, linear interpolation inbetween.
        this is the spinsight convention.

        If a pivot is provided the phase will not be changed at the pivot, but the total change accross the entire spectrum will still amount to phi1.
        """
        self.checkToPos(toPos)
        phaseValues = np.linspace(0,phi1,num=len(self.frequency))

        if noChangeOnResonance == True:
            pivot = 0
        elif pivot !=0:
            print("Using pivot for first order phase correction")
            index = self.getIndex(pivot, scale= scale)
            phaseValues = phaseValues - phaseValues[index]
            0
        if degree == True:
            phaseValues = phaseValues*np.pi/180

        self.allFid[toPos] =  [spectrum*np.exp(-1j*phaseValues) for spectrum in self.allFid[fromPos]]

def __MR_W_D_matrix(self,img,labels):
        s = sp.amax(labels)+1
        vect = self.__MR_superpixel_mean_vector(img,labels)

        adj = self.__MR_get_adj_loop(labels)

        W = sp.spatial.distance.squareform(sp.spatial.distance.pdist(vect))

        W = sp.exp(-1*W / self.weight_parameters['delta'])
        W[adj.astype(np.bool)] = 0


        D = sp.zeros((s,s)).astype(float)
        for i in range(s):
            D[i, i] = sp.sum(W[i])

        return W,D

def model(THETA, time, kplanets):
    modelo = 0.0
    if kplanets == 0:
        return 0.0
    for i in range(kplanets):
        As, P, Ac, S, C = THETA[5*i:5*(i+1)]
        A = As ** 2 + Ac ** 2
        ecc = S ** 2 + C ** 2
        w = sp.arccos(C / (ecc ** 0.5))  # longitude of periastron
        phase = sp.arccos(Ac / (A ** 0.5))
        ### test
        if S < 0:
            w = 2 * sp.pi - sp.arccos(C / (ecc ** 0.5))
        if As < 0:
            phase = 2 * sp.pi - sp.arccos(Ac / (A ** 0.5))
        ###
        per = sp.exp(P)
        freq = 2. * sp.pi / per
        M = freq * time + phase
        E = sp.array([MarkleyKESolver().getE(m, ecc) for m in M])
        f = (sp.arctan(((1. + ecc) ** 0.5 / (1. - ecc) ** 0.5) * sp.tan(E / 2.)) * 2.)
        modelo += A * (sp.cos(f + w) + ecc * sp.cos(w))
    return  modelo

def _BlackScholesCall(S, K, T, sigma, r, q):
            d1 = (log(S/K) + (r - q + (sigma**2)/2)*T)/(sigma*sqrt(T))
            d2 = d1 - sigma*sqrt(T)
            return S*exp(-q*T)*norm.cdf(d1) - K*exp(-r*T)*norm.cdf(d2)

def _BlackScholesPut(S, K, T, sigma, r, q):
            d1 = (log(S/K) + (r - q + (sigma**2)/2)*T)/(sigma*sqrt(T))
            d2 = d1 - sigma*sqrt(T)
            return  K*exp(-r*T)*norm.cdf(-d2) - S*exp(-q*T)*norm.cdf(-d1)

def _fprime(self, sigma):
        logSoverK = log(self.S/self.K)
        n12 = ((self.r + sigma**2/2)*self.T)
        numerd1 = logSoverK + n12
        d1 = numerd1/(sigma*sqrt(self.T))
        return self.S*sqrt(self.T)*norm.pdf(d1)*exp(-self.r*self.T)

def kernel(s):
    """
    Constructs a normalized discrete 1D gaussian kernel
    """
    size_grid = int(s*4)
    x= scipy.mgrid[-size_grid:size_grid+1]
    g = scipy.exp(-(x**2/float(s**2)/2.))
    return g / np.sum(g)

def alpha_m(self, V):
        """Channel gating kinetics. Functions of membrane voltage"""
        return 0.1*(V+40.0)/(1.0 - sp.exp(-(V+40.0) / 10.0))

def beta_m(self, V):
        """Channel gating kinetics. Functions of membrane voltage"""
        return 4.0*sp.exp(-(V+65.0) / 18.0)

def alpha_h(self, V):
        """Channel gating kinetics. Functions of membrane voltage"""
        return 0.07*sp.exp(-(V+65.0) / 20.0)

def beta_h(self, V):
        """Channel gating kinetics. Functions of membrane voltage"""
        return 1.0/(1.0 + sp.exp(-(V+35.0) / 10.0))

def alpha_n(self, V):
        """Channel gating kinetics. Functions of membrane voltage"""
        return 0.01*(V+55.0)/(1.0 - sp.exp(-(V+55.0) / 10.0))

def geom_std(values: t.List[float]) -> float:
    """
    Calculates the geometric standard deviation for the passed values.
    Source: https://en.wikipedia.org/wiki/Geometric_standard_deviation
    """
    import scipy.stats as stats
    import scipy as sp
    gmean = stats.gmean(values)
    return sp.exp(sp.sqrt(sp.sum([sp.log(x / gmean) ** 2 for x in values]) / len(values)))

def asian_vol_adj(atm, time2mat, tau):
    M = (2 * numpy.exp(atm * atm * time2mat) \
         - 2 * numpy.exp(atm * atm * tau) * (1.0 + atm * atm * (time2mat - tau))) / \
        ((atm ** 4) * ((time2mat - tau) ** 2))
    return numpy.sqrt(numpy.log(M) / time2mat)

def BlackSholesFormula(IsCall, S, K, Vol, Texp, Rd, Rf):
    x1 = d1(S, K, Vol, Texp, Rd, Rf )
    x2 = d2(S, K, Vol, Texp, Rd, Rf )
    y = pnorm(x1)
    res = {}

    if IsCall:
        res['Price'] = S * exp(-Rf*Texp)* x1 - K * exp(-Rd*Texp) * x2
        res['Delta'] = x1 * exp(-Rf*Texp)
    else:
        res['Price'] = K * exp(-Rd*Texp) * (1 - x2) - S * exp(-Rf*Texp) * (1 - x1)
        res['Delta'] = (x1  - 1) * exp(-Rf*Texp)
    res['Vega'] = S * sqrt(Texp) * y * exp(-Rf*Texp)
    res['Gamma'] = y * exp(-Rf*Texp)/(S*Vol* sqrt(Texp))
    return res

def KirkApprox(IsCall, F1, F2, Sigma1, Sigma2, Corr, K, Texp, r):
    FA = F1/(F2+K)
    Sigma = sqrt(Sigma1**2 + (Sigma2*F2/(F2+K))**2 - \
                       2*Corr*Sigma1*Sigma2*F2/(F2+K))
    d1 = (numpy.log(FA) + 0.5* Sigma**2 * Texp)/(Sigma*sqrt(Texp))
    d2 = d1 - Sigma*sqrt(Texp)
    x1 = scipy.stats.norm.cdf(d1)
    x2 = scipy.stats.norm.cdf(d2)
    res = {}
    if IsCall:
        res['Price'] = (F2+K)*(FA * x1 - x2) * exp(-r*Texp)
    else:
        res['Price'] = (F2+K)*((1 - x2) - FA*(1 - x1)) * exp(-r*Texp)
    return res

def BSOpt( IsCall, Spot, Strike, Vol, Texp, Rd, Rf ):
    'Standard Black-Scholes European vanilla pricing.'
    if Strike <= 1e-12 * Spot:
        if IsCall:
            return Spot * exp( -Rf * Texp )
        else:
            return 0.

    if IsCall:
        return Spot   * exp( -Rf * Texp ) * cnorm( d1( Spot, Strike, Vol, Texp, Rd, Rf ) ) \
             - Strike * exp( -Rd * Texp ) * cnorm( d2( Spot, Strike, Vol, Texp, Rd, Rf ) )
    else:
        return Strike * exp( -Rd * Texp ) * cnorm( -d2( Spot, Strike, Vol, Texp, Rd, Rf ) ) \
             - Spot   * exp( -Rf * Texp ) * cnorm( -d1( Spot, Strike, Vol, Texp, Rd, Rf ) )

def BSFwd( IsCall, Fwd, Strike, Vol, Texp, ir):
    'Standard Black-Scholes European vanilla pricing.'

    if Strike <= 1e-12 * Fwd:
        if IsCall:
            return Fwd
        else:
            return 0.
    df = exp(-ir * Texp)
    if IsCall:
        return df * (Fwd  * cnorm( fd1( Fwd, Strike, Vol, Texp ) ) \
             - Strike * cnorm( fd2( Fwd, Strike, Vol, Texp ) ))
    else:
        return df * (Strike * cnorm( -fd2( Fwd, Strike, Vol, Texp ) ) \
             - Fwd * cnorm( -fd1( Fwd, Strike, Vol, Texp ) ))

def BSFwdNormal( IsCall, Fwd, Strike, Vol, Texp, ir):
    'Standard Bachelier European vanilla pricing.'
    d = (Fwd-Strike)/Vol/sqrt(Texp)
    p = (Fwd-Strike)  * cnorm( d ) + Vol * sqrt(Texp) * pnorm(d)
    if not IsCall:
        p = p - Fwd + Strike
    return p * exp(-Texp*ir)