Predirre le vendite è il sogno più bagnato di ogni manager che si rispetti. Difatti, avere stime affidabili sull’andamento del mercato è la chiave del budgeting: più le stime sono corrette, più la programmazione diventa facile ed effettiva.
Come in ogni problema di predizione, la variabile che influisce di più sulla bontà del modello è l’informazione disponibile. Nella mia esperienza lavorativa, la quasi totalità delle aziende ha (per ora) una cultura del dato molto scarsa, se non nulla. Il data scientist si troverà dunque a dover lavorare con dati parziali, da reclamare tramite molteplici fonti, la cui affidabilità sarà semplicemente dubbia. Ovviamente, tale ostacolo è prettamente individuale e varierà da scenario a scenario: è dunque fuori dallo scopo di questo post.
Quello di cui voglio veramente parlare è come approcciare the algoritmic side of thigs . Ordunque, vi presento qui uno scoppiettante programmino per affrontare il forecasting delle serie storiche come dei veri pro.
Le basi
Il problema può essere affrontato in base a due variabili interdipendenti:
- La finestra di predizione
- Il lag delle variabili
La finestra di predizione rappresenta con quanto anticipo vogliamo predirre Y (e.g. se la predizione ha 3 mesi di anticipo, la finestra sarà di 3 mesi). Il lag delle variabili rappresenta invece la relazione temporale con cui le variabili indipendenti influenzano la dipendente. Ad esempio: le vendite di un’industria metallurgica saranno influenzate dall’andamento del PIL globale. Ma questa relazione non sarà immediata temporalmente: l’aumento di domanda del metallo sarà successivo di x mesi all’aumento del PIL.
Mentre la finestra di predizione è decisa dal Data Scientist e dal management, il lag delle variabile è fondamentalmente ignoto. Per trovarlo, in maniera ottimale, si ricorre alla Cross Validation. L’approccio che io propongo è univariato.
Restate con me, così vi spiego le funzioni per questo scopo.
L’oracolo algoritmico
Il programma è costituito da tre funzioni:
- laggedDataMat: lagga il dataset (ossia “sposta” le variabili di un determinato lag)
- bestLagRegr: trova il lag ottimale (univariato) per ogni variabile
- modelSelection: prende l’outup di bestLagRegr ed esegue una grid search per trovare l’ottimale modello predittivo
Le librerie richieste per l’esecuzione sono:
import numpy as np, pandas as pd, pickle, copy from sklearn.model_selection import cross_val_score, TimeSeriesSplit, GridSearchCV from sklearn import ensemble from sklearn.ensemble import GradientBoostingRegressor
Ecco qui di seguito le funzioni in Python:
laggedDataMat
def bestLagRegr(dataMat, maxLag, minLag, yName): # Purpose: identify the best lag for variables to be used in regression. Uses CrossValidation # (folds TimesSeriesSplit) # data - a panda DF with named columns # maxLag - the maximum lag to be tested. Integer # minLag - the minimum lag to be tested. Integer # yName - the name of the column of 'data' containing the dependent var. String data = dataMat.copy() colnames = [y for y in data.columns if y != yName] lags = range(minLag, maxLag) folds = TimeSeriesSplit(n_splits=int(np.round(data.shape[0] / 10, 0))) results = {} for col in colnames: scores = [] lags_list = [] for l in lags: varname = col + str('_lag') + str(l) data[varname] = data[col].shift(l) YX = data[[yName, varname]] YX = YX.dropna().as_matrix() # Build regressor and estimate metric of prediction performance with CV regr = ensemble.GradientBoostingRegressor(learning_rate=0.01, max_depth=1, n_estimators=500) perform = cross_val_score(regr, X=YX[:, 1].reshape(-1, 1), y=YX[:, 0], cv=folds, scoring='neg_median_absolute_error') # Store scores result and lags scores.append(np.median(perform)) lags_list.append(l) # Calculate best score, corresponding best lag and store it in a dictionary, containing all colnames best_score = max(scores) best_lag = lags_list[scores.index(best_score)] results[col] = [best_lag, best_score] return(results)
bestLagRegr
def laggedDataMat(dataMat, yName, lagsDict): # Purpose: build a lagged DF # dataMat: the unlagged DF with named columns # yName: name of the dependent var. String # lagsDict: dictionary produced with 'bestLagRegr, containing: # - keys: column names of dataMat # - elements: lists with lag, CV score # # Output: a panda DataFrame with columns order sorted alphabetically # Initialize empty DF df = pd.DataFrame(index=dataMat.index) # Set dependent var df[yName] = dataMat[[yName]] # Creating and adding the lagged vars for colName in lagsDict.keys(): l = lagsDict[colName][0] colNameLag = colName + str('_lag') + str(l) df[colNameLag] = dataMat[[colName]].shift(l) df = df.sort_index(axis=1) return(df)
modelSelection
def modelSelection(maxLag, data, depVar, toSave, pSpace = None, alpha = 0.95): if pSpace is None: pSpace = dict(n_estimators=list(range(5, 2000, 10)), learning_rate=list(np.arange(0.001, 1, 0.1)), max_depth=list(range(1, 3, 1))) lags = range(1, maxLag) results = dict() for lagMin in lags: print('Esimating model for lag: ', lagMin) lagAnalysis = bestLagRegr(data, maxLag, lagMin, depVar) lagMat = laggedDataMat(data, depVar, lagAnalysis) lagMat = lagMat.dropna() trainY = np.ravel(lagMat[depVar]) lagMat = lagMat.drop([depVar], 1) folds = TimeSeriesSplit(n_splits=int(round(lagMat.shape[0] / 10, 0))) model = GradientBoostingRegressor(loss='ls') regr = GridSearchCV(estimator=model, param_grid=pSpace, scoring='neg_mean_squared_error', cv=folds) regr.fit(lagMat, trainY) modelName = toSave + '/' + 'modelOI' + '_lag' + str(lagMin) + '.sav' pickle.dump(regr.best_estimator_, open(modelName, 'wb')) temp = dict(BestModel=regr.best_estimator_, Score=regr.best_score_, Lags=lagAnalysis) regrQUpper = copy.deepcopy(regr.best_estimator_) regrQUpper.set_params(loss='quantile', alpha=alpha) regrQUpper.fit(lagMat, trainY) temp['QUpper'] = regrQUpper regrQLower = copy.deepcopy(regr.best_estimator_) regrQLower.set_params(loss='quantile', alpha=(1-alpha)) regrQLower.fit(lagMat, trainY) temp['QLower'] = regrQLower key = 'Lag ' + str(lagMin) results[key] = temp return(results)