Erstelle mehrere Spalten in Pandas Dataframe von einer Funktion

Question

Ich bin ein Python-Neuling, also hoffe ich, dass meine zwei Fragen klar und vollständig sind. Ich habe den tatsächlichen Code und einen Testdatensatz im CSV-Format unten veröffentlicht.

Ich konnte den folgenden Code (meist mit Hilfe der StackOverflow-Mitwirkenden) erstellen, um die implizite Volatilität eines Optionskontrakts mit der Newton-Raphson-Methode zu berechnen. Der Prozess berechnet Vega bei der Bestimmung der impliziten Volatilität. Obwohl ich eine neue DataFrame-Spalte für Implied Volatility mithilfe der Pandas DataFrame-Methode apply erstellen kann, kann ich keine zweite Spalte für Vega erstellen. Gibt es eine Möglichkeit, zwei separate DataFrame-Spalten zu erstellen, wenn die Funktion IV & Vega zusammen zurückgibt?

Ich habe versucht:

• return iv, vega von Funktion
• df[['myIV', 'Vega']] = df.apply(newtonRap, axis=1)
• Got ValueError: Shape of passed values is (56, 2), indices imply (56, 13)

auch versucht:

• return iv, vega von Funktion
• df['myIV'], df['Vega'] = df.apply(newtonRap, axis=1)
• Got ValueError: Shape of passed values is (56, 2), indices imply (56, 13)

Zusätzlich ist der Berechnungsvorgang langsam. Ich habe numba importiert und den @ jit (nogil = True) Decorator implementiert, aber ich sehe nur eine Leistungsverbesserung von 25%. Der Testdatensatz ist der Leistungstest hat fast 900.000 Datensätze. Die Laufzeit beträgt 2 Stunden und 9 Minuten ohne Numba oder mit Numba, aber ohne Nogil = True. Die Laufzeit bei Verwendung von numba und @jit (nogil = True) beträgt 1 Stunde und 32 Minuten. Kann ich es besser machen?

from datetime import datetime 
from math import sqrt, pi, log, exp, isnan 
from scipy.stats import norm 
from numba import jit 


# dff = Daily Fed Funds (Posted rate is usually one day behind) 
dff = pd.read_csv('https://research.stlouisfed.org/fred2/data/DFF.csv', parse_dates=[0], index_col='DATE') 
rf = float('%.4f' % (dff['VALUE'][-1:][0]/100)) 
# rf = .0015      # Get Fed Funds Rate https://research.stlouisfed.org/fred2/data/DFF.csv 
tradingMinutesDay = 450    # 7.5 hours per day * 60 minutes per hour 
tradingMinutesAnnum = 113400  # trading minutes per day * 252 trading days per year 
cal = USFederalHolidayCalendar() # Load US Federal holiday calendar 


@jit(nogil=True)        # nogil=True arg improves performance by 25% 
def newtonRap(row): 
    """Estimate Implied Volatility (IV) using Newton-Raphson method 

    :param row (dataframe): Options contract params for function 
     TimeStamp (datetime): Close date 
     Expiry (datetime): Option contract expiration date 
     Strike (float): Option strike 
     OptType (object): 'C' for call; 'P' for put 
     RootPrice (float): Underlying close price 
     Bid (float): Option contact closing bid 
     Ask (float): Option contact closing ask 

    :return: 
     float: Estimated implied volatility 
    """ 
    if row['Bid'] == 0.0 or row['Ask'] == 0.0 or row['RootPrice'] == 0.0 or row['Strike'] == 0.0 or \ 
     row['TimeStamp'] == row['Expiry']: 
     iv, vega = 0.0, 0.0   # Set iv and vega to zero if option contract is invalid or expired 
    else: 
     # dte (Days to expiration) uses pandas bdate_range method to determine the number of business days to expiration 
     # minus USFederalHolidays minus constant of 1 for the TimeStamp date 
     dte = float(len(pd.bdate_range(row['TimeStamp'], row['Expiry'])) - 
        len(cal.holidays(row['TimeStamp'], row['Expiry']).to_pydatetime()) - 1) 
     mark = (row['Bid'] + row['Ask'])/2 
     cp = 1 if row['OptType'] == 'C' else -1 
     S = row['RootPrice'] 
     K = row['Strike'] 
     # T = the number of trading minutes to expiration divided by the number of trading minutes in year 
     T = (dte * tradingMinutesDay)/tradingMinutesAnnum 
     # TODO get dividend value 
     d = 0.00 
     iv = sqrt(2 * pi/T) * mark/S  # Closed form estimate of IV Brenner and Subrahmanyam (1988) 
     vega = 0.0 
     for i in range(1, 100): 
      d1 = (log(S/K) + T * (rf - d + iv ** 2/2))/(iv * sqrt(T)) 
      d2 = d1 - iv * sqrt(T) 
      vega = S * norm.pdf(d1) * sqrt(T) 
      model = cp * S * norm.cdf(cp * d1) - cp * K * exp(-rf * T) * norm.cdf(cp * d2) 
      iv -= (model - mark)/vega 
      if abs(model - mark) < 1.0e-9: 
       break 
     if isnan(iv) or isnan(vega): 
      iv, vega = 0.0, 0.0 
    # TODO Return vega with iv if add'l pandas column possible 
    # return iv, vega 
    return iv 


if __name__ == "__main__": 
    # test function from baseline data 
    get_csv = True 

    if get_csv: 
     csvHeaderList = ['TimeStamp', 'OpraSymbol', 'RootSymbol', 'Expiry', 'Strike', 'OptType', 'RootPrice', 'Last', 
         'Bid', 'Ask', 'Volume', 'OpenInt', 'IV'] 
     fileName = 'C:/tmp/test-20150930-56records.csv' 
     df = pd.read_csv(fileName, parse_dates=[0, 3], names=csvHeaderList) 
    else: 
     pass 

    start = datetime.now() 
    # TODO Create add'l pandas dataframe column, if possible, for vega 
    # df[['myIV', 'Vega']] = df.apply(newtonRap, axis=1) 
    # df['myIV'], df['Vega'] = df.apply(newtonRap, axis=1) 
    df['myIV'] = df.apply(newtonRap, axis=1) 
    end = datetime.now() 
    print end - start 

Testdaten: C: /tmp/test-20150930-56records.csv

2015.09.30 16: 00: 00, AAPL151016C00109000, AAPL, 2015.10.16 16.00 Uhr: 00,109, C, 109,95,3,46,3,6,3,7,1565,1290,0,3497 2015-09-30 16: 00: 00, AAPL151016P00109000, AAPL, 2015-10-16 16: 00: 00,109, P, 109,95,2,4, 2.34,2,42,3790,308,0,3,06,10 2015-09-30 16: 00: 00, AAPL151016C00110000, AAPL, 2015-10-16 16: 00: 00,110, C, 109,95,3,2,86,3,10217,28850, 0.3288 2015-09-30 16: 00: 00, AAPL151016P00110000, AAPL, 2015-10-16 16: 00: 00,110, P, 109,95,2,81,2,74,2,8,12113,44427,0,3029 2015-09-30 16 : 00: 00, AAPL151016C00111000, AAPL, 2015-10-16 1 6: 00: 00,111, C, 109,95,2,35,2,44,2,45,6674,2318,0,3187 2015-09-30 16: 00: 00, AAPL151016P00111000, AAPL, 2015-10-16 16: 00: 00,111, P, 109,95,3,2,3,1,3,25,2031,3773,0,2926 2015-09-30 16: 00: 00, AAPL151120C00110000, AAPL, 2015-11-20 16: 00: 00,110, C, 109,95,5,9,5,7,5,95, 5330,17112,0,3635 2015-09-30 16: 00: 00, AAPL151120P00110000, AAPL, 2015-11-20 16: 00: 00,110, P, 109,95,6,15,6,1,6,3,3724,15704,0,3842

Answer 1

Wenn ich Sie richtig verstehe, sollten Sie eine Serie von Ihrer Funktion zurückgeben.Etwas wie:

return pandas.Series({"IV": iv, "Vega": vega}) 

Wenn Sie das Ergebnis in neue Spalten des gleichen Eingangsdatenrahmen setzen wollen, dann tun Sie nur:

df[["IV", "Vega"]] = df.apply(newtonRap, axis=1) 

Answer 2

Was die Leistung mit numba betrifft, numba doesn‘ Ich weiß nichts über Pandas-Datenframes und kann Operationen auf ihnen nicht zu schnellem Maschinencode kompilieren. Am besten profilieren Sie sich, welcher Teil Ihrer Methode langsam ist (z. B. line_profiler), und laden Sie dann diesen Teil auf eine andere Methode, die die Eingaben unter Verwendung der .values-Attribute der Datenframe-Spalten erstellt. Dadurch erhalten Sie Zugriff auf die zugrundeliegende Zahl Array. Sonst wird numba hauptsächlich im "Objektmodus" betrieben (siehe numba glossary) und wird die Leistung nicht drastisch verbessern

Answer 3

Der Trick, Code zu vektorisieren, besteht darin, nicht in Zeilen zu denken, sondern stattdessen in Spalten zu denken .

Ich habe fast diese Arbeit (ich werde versuchen, es später zu beenden), aber Sie wollen etwas entlang der Linien von, dies zu tun:

from datetime import datetime 
from math import sqrt, pi, log, exp, isnan 
from numpy import inf, nan 
from scipy.stats import norm 
import pandas as pd 
from pandas import Timestamp 
from pandas.tseries.holiday import USFederalHolidayCalendar 

# Initial parameters 
rf = .0015       # Get Fed Funds Rate https://research.stlouisfed.org/fred2/data/DFF.csv 
tradingMinutesDay = 450    # 7.5 hours per day * 60 minutes per hour 
tradingMinutesAnnum = 113400  # trading minutes per day * 252 trading days per year 
cal = USFederalHolidayCalendar() # Load US Federal holiday calendar 
two_pi = 2 * pi      # 2 * Pi (to reduce computations) 
threshold = 1.0e-9     # convergence threshold. 

# Create sample data: 
col_order = ['TimeStamp', 'OpraSymbol', 'RootSymbol', 'Expiry', 'Strike', 'OptType', 'RootPrice', 'Last', 'Bid', 'Ask', 'Volume', 'OpenInt', 'IV'] 
df = pd.DataFrame({'Ask': {0: 3.7000000000000002, 1: 2.4199999999999999, 2: 3.0, 3: 2.7999999999999998, 4: 2.4500000000000002, 5: 3.25, 6: 5.9500000000000002, 7: 6.2999999999999998}, 
        'Bid': {0: 3.6000000000000001, 1: 2.3399999999999999, 2: 2.8599999999999999, 3: 2.7400000000000002, 4: 2.4399999999999999, 5: 3.1000000000000001, 6: 5.7000000000000002, 7: 6.0999999999999996}, 
        'Expiry': {0: Timestamp('2015-10-16 16:00:00'), 1: Timestamp('2015-10-16 16:00:00'), 2: Timestamp('2015-10-16 16:00:00'), 3: Timestamp('2015-10-16 16:00:00'), 4: Timestamp('2015-10-16 16:00:00'), 5: Timestamp('2015-10-16 16:00:00'), 6: Timestamp('2015-11-20 16:00:00'), 7: Timestamp('2015-11-20 16:00:00')}, 
        'IV': {0: 0.3497, 1: 0.3146, 2: 0.3288, 3: 0.3029, 4: 0.3187, 5: 0.2926, 6: 0.3635, 7: 0.3842}, 
        'Last': {0: 3.46, 1: 2.34, 2: 3.0, 3: 2.81, 4: 2.35, 5: 3.20, 6: 5.90, 7: 6.15}, 
        'OpenInt': {0: 1290.0, 1: 3087.0, 2: 28850.0, 3: 44427.0, 4: 2318.0, 5: 3773.0, 6: 17112.0, 7: 15704.0}, 
        'OpraSymbol': {0: 'AAPL151016C00109000', 1: 'AAPL151016P00109000', 2: 'AAPL151016C00110000', 3: 'AAPL151016P00110000', 4: 'AAPL151016C00111000', 5: 'AAPL151016P00111000', 6: 'AAPL151120C00110000', 7: 'AAPL151120P00110000'}, 
        'OptType': {0: 'C', 1: 'P', 2: 'C', 3: 'P', 4: 'C', 5: 'P', 6: 'C', 7: 'P'}, 
        'RootPrice': {0: 109.95, 1: 109.95, 2: 109.95, 3: 109.95, 4: 109.95, 5: 109.95, 6: 109.95, 7: 109.95}, 
        'RootSymbol': {0: 'AAPL', 1: 'AAPL', 2: 'AAPL', 3: 'AAPL', 4: 'AAPL', 5: 'AAPL', 6: 'AAPL', 7: 'AAPL'}, 
        'Strike': {0: 109.0, 1: 109.0, 2: 110.0, 3: 110.0, 4: 111.0, 5: 111.0, 6: 110.0, 7: 110.0}, 
        'TimeStamp': {0: Timestamp('2015-09-30 16:00:00'), 1: Timestamp('2015-09-30 16:00:00'), 2: Timestamp('2015-09-30 16:00:00'), 3: Timestamp('2015-09-30 16:00:00'), 4: Timestamp('2015-09-30 16:00:00'), 5: Timestamp('2015-09-30 16:00:00'), 6: Timestamp('2015-09-30 16:00:00'), 7: Timestamp('2015-09-30 16:00:00')}, 
        'Volume': {0: 1565.0, 1: 3790.0, 2: 10217.0, 3: 12113.0, 4: 6674.0, 5: 2031.0, 6: 5330.0, 7: 3724.0}}) 
df = df[col_order] 

# Vectorize columns 
df['mark'] = (df.Bid + df.Ask)/2 
df['cp'] = df.OptType.map({'C': 1, 'P': -1}) 
df['Log_S_K'] = (df.RootPrice/df.Strike).apply(log) 
df['divs'] = 0 # TODO: Get dividend value. 
df['vega'] = 0. 
df['converged'] = False 

# Vectorized datetime calculations 
date_pairs = set(zip(df.TimeStamp, df.Expiry)) 
total_days = {(t1, t2): len(pd.bdate_range(t1, t2)) 
         for t1, t2 in date_pairs} 
hols = {(t1, t2): len(cal.holidays(t1, t2).to_pydatetime()) 
        for t1, t2 in date_pairs} 
del date_pairs 

df['total_days'] = [total_days.get((t1, t2)) 
        for t1, t2 in zip(df.TimeStamp, df.Expiry)] 
df['hols'] = [hols.get((t1, t2)) 
       for t1, t2 in zip(df.TimeStamp, df.Expiry)] 
df['days_to_exp'] = df.total_days - df.hols - 1 
df.loc[df.days_to_exp < 0, 'days_to_exp'] = 0 # Min zero. 
df.drop(['total_days', 'hols'], axis='columns', inplace=True) 
df['years_to_expiry'] = (df.days_to_exp * tradingMinutesDay/tradingMinutesAnnum) 

# Initial implied vol 'guess' 
df['implied_vol'] = (two_pi/df.years_to_expiry) ** 0.5 * df.mark/df.RootPrice 

for i in xrange(100): # range(100) in Python 3.x 
    # Create mask of options where the vol has not converged. 
    mask = [not c for c in df.converged.values] 
    if df.converged.all(): 
     break 

    # Aliases. 
    data = df.loc[mask, :] 
    cp = data.cp 
    mark = data.mark 
    S = data.RootPrice 
    K = data.Strike 
    d = data.divs 
    T = data.years_to_expiry 
    log_S_K = data.Log_S_K 
    iv = data.implied_vol 

    # Calcs. 
    d1 = (log_S_K + T * (rf - d + .5 * iv ** 2))/(iv * T ** 0.5) 
    d2 = d1 - iv * T ** 0.5 
    df.loc[mask, 'vega'] = vega = S * d1.apply(norm.pdf) * T ** 0.5 
    model = cp * (S * (cp * d1).apply(norm.cdf) 
        - K * (-rf * T).apply(exp) * (cp * d2).apply(norm.cdf)) 
    iv_delta = (model - mark)/vega 
    df.loc[mask, 'implied_vol'] = iv - iv_delta 

    # Clean-up and check for convergence. 
    df.loc[df.implied_vol < 0, 'implied_vol'] = 0 
    idx = model[(model - mark).abs() < threshold].index 
    df.ix[idx, 'converged'] = True 
    df.loc[:, 'implied_vol'].fillna(0, inplace=True) 
    df.loc[:, 'implied_vol'].replace([inf, -inf], nan, inplace=True) 
    df.loc[:, 'vega'].fillna(0, inplace=True) 
    df.loc[:, 'vega'].replace([inf, -inf], nan, inplace=True)