1 '''
2 (c) 2011, 2012 Georgia Tech Research Corporation
3 This source code is released under the New BSD license. Please see
4 http://wiki.quantsoftware.org/index.php?title=QSTK_License
5 for license details.
6
7 Created on September, 12, 2011
8
9 @author: Who?
10 @contact:
11 @summary: Time Series utilities.
12 '''
13
14 import cPickle
15 import math
16 import datetime as dt
17 import numpy as np
18 from pylab import *
19 from pandas import *
20 from qstkutil import dateutil
21 from math import sqrt, log
22 from copy import deepcopy
23
24 import random as rand
25
26 from qstkutil import DataAccess as da
27 from qstkutil import dateutil as du
28
30
31
32 nd=deepcopy(funds)
33 nd[0]=0
34
35 for i in range(1,len(funds)):
36 nd[i]=funds[i]/funds[i-1]-1
37 return(nd)
38
47
49 rets=daily(funds)
50 x=0
51 years=dateutil.getYears(funds)
52 averages=[]
53 for year in years:
54 months=dateutil.getMonths(funds,year)
55 for month in months:
56 avg=0
57 count=0
58 days=dateutil.getDays(funds,year,month)
59 for day in days:
60 avg+=rets[x]
61 x+=1
62 count+=1
63 averages.append(float(avg)/count)
64 return(averages)
65
67 """
68 @summary Removes NaNs from a 2D array by scanning forward in the
69 1st dimension. If a cell is NaN, the value above it is carried forward.
70 @param nd: the array to fill forward
71 @return the array is revised in place
72 """
73 for col in range(nd.shape[1]):
74 for row in range(1,nd.shape[0]):
75 if math.isnan(nd[row,col]):
76 nd[row,col] = nd[row-1,col]
77
79 """
80 @summary Removes NaNs from a 2D array by scanning backward in the
81 1st dimension. If a cell is NaN, the value above it is carried backward.
82 @param nd: the array to fill backward
83 @return the array is revised in place
84 """
85 for col in range(nd.shape[1]):
86 for row in range(nd.shape[0]-2,-1,-1):
87 if math.isnan(nd[row,col]):
88 nd[row,col] = nd[row+1,col]
89
91 """
92 @summary Computes stepwise (usually daily) returns relative to 0, where
93 0 implies no change in value.
94 @return the array is revised in place
95 """
96 nd[1:,:] = (nd[1:,:]/nd[0:-1]) - 1
97 nd[0,:] = np.zeros(nd.shape[1])
98
100 """
101 @summary Computes stepwise (usually daily) returns relative to 1, where
102 1 implies no change in value.
103 @param nd: the array to fill backward
104 @return the array is revised in place
105 """
106 nd[1:,:] = (nd[1:,:]/nd[0:-1])
107 nd[0,:] = np.ones(nd.shape[1])
108
110 """
111 @summary Computes stepwise (usually daily) returns relative to 1, where
112 1 implies no change in value.
113 @param nd: the array to fill backward
114 @return the array is revised in place
115 """
116
117 nd[0,:] = 100
118 for i in range(1,nd.shape[0]):
119 nd[i,:] = nd[i-1,:] * nd[i,:]
120
121
123 """
124 @summary Computes stepwise (usually daily) logarithmic returns.
125 @param nd: the array to fill backward
126 @return the array is revised in place
127 """
128 returnize1(nd)
129 nd = np.log(nd)
130 return nd
131
133 d=daily(funds)
134 avg=float(sum(d))/len(d)
135 std=0
136 for a in d:
137 std=std+float((float(a-avg))**2)
138 std=sqrt(float(std)/(len(d)-1))
139 return(avg/std)
140
142 funds2=[]
143 for date in funds.index:
144 if(date.year==year):
145 funds2.append(funds[date])
146 return(getRatio(funds2))
147
149 """
150 @summary Returns the daily Sharpe ratio of the returns.
151 @param naRets: 1d numpy array or list of daily returns
152 @param fFreeReturn: risk free returns, default is 3%
153 @return Annualized rate of return, not converted to percent
154 """
155
156 fDev = np.std( naRets - 1, axis=0 )
157 fMean = np.mean( naRets - 1, axis=0 )
158
159 ''' Convert to yearly standard deviation '''
160 fSharpe = (fMean * 252 - fFreeReturn) / ( fDev * sqrt(252) )
161
162
163
164 return fSharpe
165
167 """
168 @summary Returns the rate of return annualized. Assumes len(naRets) is number of days.
169 @param naRets: 1d numpy array or list of daily returns
170 @return Annualized rate of return, not converted to percent
171 """
172
173 ''' Calculate final value of investment of 1.0 '''
174 fInv = 1.0
175 for fReturn in naRets:
176 fInv = fInv * fReturn
177
178 fRorYtd = fInv - 1.0
179
180 print ' RorYTD =', fInv, 'Over days:', len(naRets)
181
182 return ( (1.0 + fRorYtd)**( 1.0/(len(naRets)/365.0) ) ) - 1.0
183
185 """
186 @summary Reindexes a DataMatrix price array and returns the new periodic returns.
187 @param dmPrice: DataMatrix of stock prices
188 @param sOffset: Offset string to use, choose from _offsetMap in pandas/core/datetools.py
189 e.g. 'EOM', 'WEEKDAY', 'W@FRI', 'A@JAN'. Or use a pandas DateOffset.
190 """
191
192 ''' Could possibly use DataMatrix.asfreq here '''
193 ''' Use pandas DateRange to create the dates we want, use 4:00 '''
194 drNewRange = DateRange(dmPrice.index[0], dmPrice.index[-1], timeRule=sOffset) + DateOffset(hours=16)
195
196 dmPrice = dmPrice.reindex( drNewRange, method='ffill' )
197
198 returnize1( dmPrice.values )
199
200 ''' Do not leave return of 1.0 for first time period: not accurate '''
201 return dmPrice[1:]
202
204 """
205 @summary Reindexes returns using the cumulative product. E.g. if returns are 1.5 and 1.5, a period of 2 will
206 produce a 2-day return of 2.25. Note, these must be returns centered around 1.
207 @param naRets: Daily returns of the various stocks (using returnize1)
208 @param lPeriod: New target period.
209 @note: Note that this function does not track actual weeks or months, it only approximates with trading days.
210 You can use 5 for week, or 21 for month, etc.
211 """
212 naCumData = np.cumprod(naRets, axis=0)
213
214 lNewRows =(naRets.shape[0]-1) / (lPeriod)
215 ''' We compress data into height / lPeriod + 1 new rows '''
216 for i in range( lNewRows ):
217 lCurInd = -1 - i*lPeriod
218 ''' Just hold new data in same array, new return is cumprod on day x / cumprod on day x-lPeriod '''
219 naCumData[-1 - i,:] = naCumData[lCurInd,:] / naCumData[lCurInd - lPeriod,:]
220 ''' Select new returns from end of cumulative array '''
221
222 return naCumData[-lNewRows:, ]
223
224
225 -def getOptPort( naRets, fTarget, lPeriod=1, naLower=None, naUpper=None, lNagDebug=0 ):
226 """
227 @summary Returns the Markowitz optimum portfolio for a specific return.
228 @param naRets: Daily returns of the various stocks (using returnize1)
229 @param fTarget: Target return, i.e. 0.04 = 4% per period
230 @param lPeriod: Period to compress the returns to, e.g. 7 = weekly
231 @param naLower: List of floats which corresponds to lower portfolio% for each stock
232 @param naUpper: List of floats which corresponds to upper portfolio% for each stock
233 @return tuple: (weights of portfolio, min possible return, max possible return)
234 """
235
236 ''' Attempt to import library '''
237 try:
238 pass
239 import nagint as nag
240 except ImportError:
241 print 'Could not import NAG library, make sure nagint.so is in your python path'
242 return ([], 0, 0)
243
244 ''' Get number of stocks '''
245 lStocks = naRets.shape[1]
246
247 ''' If period != 1 we need to restructure the data '''
248 if( lPeriod != 1 ):
249 naRets = getReindexedRets( naRets, lPeriod)
250
251 ''' Calculate means and covariance '''
252 naAvgRets = np.average( naRets, axis=0 )
253 naCov = np.cov( naRets, rowvar=False )
254
255 ''' Special case for None == fTarget, simply return average returns and cov '''
256 if( fTarget is None ):
257 return naAvgRets, np.std(naRets, axis=0)
258
259 ''' Calculate upper and lower limits of variables as well as constraints '''
260 if( naUpper is None ):
261 naUpper = np.ones( lStocks )
262
263 if( naLower is None ):
264 naLower = np.zeros( lStocks )
265 ''' Two extra constraints for linear conditions, result = desired return, and sum of weights = 1 '''
266 naUpper = np.append( naUpper, [fTarget, 1.0] )
267 naLower = np.append( naLower, [fTarget, 1.0] )
268
269 ''' Initial estimate of portfolio '''
270 naInitial = np.array([1.0/lStocks]*lStocks)
271
272 ''' Set up constraints matrix, composed of expected returns in row one, unity row in row two '''
273 naConstraints = np.vstack( (naAvgRets, np.ones(lStocks)) )
274
275 ''' Get portfolio weights, last entry in array is actually variance '''
276 try:
277 naReturn = nag.optPort( naConstraints, naLower, naUpper, naCov, naInitial, lNagDebug )
278 except RuntimeError:
279 print 'NAG Runtime error with target: %.02lf'%(fTarget)
280 return ( naInitial, sqrt( naCov[0][0] ) )
281
282 ''' Calculate stdev of entire portfolio to return, what NAG returns is slightly different '''
283 fPortDev = np.std( np.dot(naRets, naReturn[0,0:-1]) )
284
285 ''' Show difference between above stdev and sqrt NAG covariance, possibly not taking correlation into account '''
286
287
288 ''' Return weights and stdDev of portfolio. note again the last value of naReturn is NAG's reported variance '''
289 return (naReturn[0,0:-1], fPortDev)
290
291
292
294 """
295 @summary Returns the range of possible returns with upper and lower bounds on the portfolio participation
296 @param naRets: Expected returns
297 @param naLower: List of lower percentages by stock
298 @param naUpper: List of upper percentages by stock
299 @return tuple containing (fMin, fMax)
300 """
301
302 ''' Calculate theoretical minimum and maximum theoretical returns '''
303 fMin = 0
304 fMax = 0
305
306 naAvgRets = np.average( naRets, axis=0 )
307 naSortInd = naAvgRets.argsort()
308
309 ''' First add the lower bounds on portfolio participation '''
310 for i, fRet in enumerate(naAvgRets):
311 fMin = fMin + fRet*naLower[i]
312 fMax = fMax + fRet*naLower[i]
313
314
315 ''' Now calculate minimum returns, allocate the max possible in worst performing equities '''
316 ''' Subtract min since we have already counted it '''
317 naUpperAdd = naUpper - naLower
318 fTotalPercent = np.sum(naLower[:])
319 for i, lInd in enumerate(naSortInd):
320 fRetAdd = naUpperAdd[lInd] * naAvgRets[lInd]
321 fTotalPercent = fTotalPercent + naUpperAdd[lInd]
322 fMin = fMin + fRetAdd
323
324 ''' Check if this additional percent puts us over the limit '''
325 if fTotalPercent > 1.0:
326 fMin = fMin - naAvgRets[lInd] * (fTotalPercent - 1.0)
327 break
328
329 ''' Repeat for max, just reverse the sort, i.e. high to low '''
330 naUpperAdd = naUpper - naLower
331 fTotalPercent = np.sum(naLower[:])
332 for i, lInd in enumerate(naSortInd[::-1]):
333 fRetAdd = naUpperAdd[lInd] * naAvgRets[lInd]
334 fTotalPercent = fTotalPercent + naUpperAdd[lInd]
335 fMax = fMax + fRetAdd
336
337 ''' Check if this additional percent puts us over the limit '''
338 if fTotalPercent > 1.0:
339 fMax = fMax - naAvgRets[lInd] * (fTotalPercent - 1.0)
340 break
341
342 return (fMin, fMax)
343
344
345 -def getFrontier( naRets, lRes=100, fUpper=0.2, fLower=0.00):
346 """
347 @summary Generates an efficient frontier based on average returns.
348 @param naRets: Array of returns to use
349 @param lRes: Resolution of the curve, default=100
350 @param fUpper: Upper bound on portfolio percentage
351 @param fLower: Lower bound on portfolio percentage
352 @return tuple containing (lfReturn, lfStd, lnaPortfolios)
353 lfReturn: List of returns provided by each point
354 lfStd: list of standard deviations provided by each point
355 lnaPortfolios: list of numpy arrays containing weights for each portfolio
356 """
357
358 ''' Limit/enforce percent participation '''
359 naUpper = np.ones(naRets.shape[1]) * fUpper
360 naLower = np.ones(naRets.shape[1]) * fLower
361
362 (fMin, fMax) = getRetRange( naRets, naLower, naUpper )
363
364 ''' Try to avoid intractible endpoints due to rounding errors '''
365 fMin *= 1.0000001
366 fMax *= 0.9999999
367
368 ''' Calculate target returns from min and max '''
369 lfReturn = []
370 for i in range(lRes):
371 lfReturn.append( (fMax - fMin) * i / (lRes - 1) + fMin )
372
373
374 lfStd = []
375 lnaPortfolios = []
376
377 ''' Call the function lRes times for the given range, use 1 for period '''
378 for fTarget in lfReturn:
379 (naWeights, fStd) = getOptPort( naRets, fTarget, 1, naUpper=naUpper, naLower=naLower )
380 lfStd.append(fStd)
381 lnaPortfolios.append( naWeights )
382
383 ''' plot frontier '''
384 '''plt.plot( lfStd, lfReturn )
385 plt.plot( np.std( naRets, axis=0 ), np.average( naRets, axis=0 ), 'g+', markersize=10 )
386 #plt.show()'''
387
388 return (lfReturn, lfStd, lnaPortfolios)
389
390
391 -def stockFilter( dmPrice, dmVolume, fNonNan=0.95, fPriceVolume=100*1000 ):
392 """
393 @summary Returns the list of stocks filtered based on various criteria.
394 @param dmPrice: DataMatrix of stock prices
395 @param dmVolume: DataMatrix of stock volumes
396 @param fNonNan: Optional non-nan percent, default is .95
397 @param fPriceVolume: Optional price*volume, default is 100,000
398 @return list of stocks which meet the criteria
399 """
400
401 lsRetStocks = list( dmPrice.columns )
402
403 for sStock in dmPrice.columns:
404 fValid = 0.0
405 print sStock
406 ''' loop through all dates '''
407 for dtDate in dmPrice.index:
408 ''' Count null (nan/inf/etc) values '''
409 fPrice = dmPrice[sStock][dtDate]
410 if( not isnull(fPrice) ):
411 fValid = fValid + 1
412 ''' else test price volume '''
413 fVol = dmVolume[sStock][dtDate]
414 if( not isnull(fVol) and fVol * fPrice < fPriceVolume ):
415 lsRetStocks.remove( sStock )
416 break
417
418 ''' Remove if too many nan values '''
419 if( fValid / len(dmPrice.index) < fNonNan and sStock in lsRetStocks ):
420 lsRetStocks.remove( sStock )
421
422 return lsRetStocks
423
424
425 -def getRandPort( lNum, dtStart=None, dtEnd=None, lsStocks=None, dmPrice=None, dmVolume=None, bFilter=True, fNonNan=0.95, fPriceVolume=100*1000, lSeed=None ):
426 """
427 @summary Returns a random portfolio based on certain criteria.
428 @param lNum: Number of stocks to be included
429 @param dtStart: Start date for portfolio
430 @param dtEnd: End date for portfolio
431 @param lsStocks: Optional list of ticker symbols, if not provided all symbols will be used
432 @param bFilter: If False, stocks are not filtered by price or volume data, simply return random Portfolio.
433 @param dmPrice: Optional price data, if not provided, data access will be queried
434 @param dmVolume: Optional volume data, if not provided, data access will be queried
435 @param fNonNan: Optional non-nan percent for filter, default is .95
436 @param fPriceVolume: Optional price*volume for filter, default is 100,000
437 @warning: Does not work for all sets of optional inputs, e.g. if you don't include dtStart, dtEnd, you need
438 to include dmPrice/dmVolume
439 @return list of stocks which meet the criteria
440 """
441
442 if( lsStocks is None ):
443 if( dmPrice is None and dmVolume is None ):
444 norObj = da.DataAccess('Norgate')
445 lsStocks = norObj.get_all_symbols()
446 elif( not dmPrice is None ):
447 lsStocks = list(dmPrice.columns)
448 else:
449 lsStocks = list(dmVolume.columns)
450
451 if( dmPrice is None and dmVolume is None and bFilter == True ):
452 norObj = da.DataAccess('Norgate')
453 ldtTimestamps = du.getNYSEdays( dtStart, dtEnd, dt.timedelta(hours=16) )
454
455 ''' if dmPrice and dmVol are provided then we don't query it every time '''
456 bPullPrice = False
457 bPullVol = False
458 if( dmPrice is None ):
459 bPullPrice = True
460 if( dmVolume is None ):
461 bPullVol = True
462
463 ''' Default seed (none) uses system clock '''
464 rand.seed(lSeed)
465 lsRetStocks = []
466
467 ''' Loop until we have enough randomly selected stocks '''
468 llRemainingIndexes = range(0,len(lsStocks))
469 lsValid = None
470 while( len(lsRetStocks) != lNum ):
471
472 lsCheckStocks = []
473 for i in range( lNum - len(lsRetStocks) ):
474 lRemaining = len(llRemainingIndexes)
475 if( lRemaining == 0 ):
476 print 'Error in getRandPort: ran out of stocks'
477 return lsRetStocks
478
479 ''' Pick a stock and remove it from the list of remaining stocks '''
480 lPicked = rand.randint(0, lRemaining-1)
481 lsCheckStocks.append( lsStocks[ llRemainingIndexes.pop(lPicked) ] )
482
483 ''' If bFilter is false, simply return our first list of stocks, don't check prive/vol '''
484 if( not bFilter ):
485 return sorted(lsCheckStocks)
486
487
488 ''' Get data if needed '''
489 if( bPullPrice ):
490 dmPrice = norObj.get_data( ldtTimestamps, lsCheckStocks, 'close' )
491
492 ''' Get data if needed '''
493 if( bPullVol ):
494 dmVolume = norObj.get_data( ldtTimestamps, lsCheckStocks, 'volume' )
495
496 ''' Only query this once if data is provided, else query every time with new data '''
497 if( lsValid is None or bPullVol or bPullPrice ):
498 lsValid = stockFilter(dmPrice, dmVolume, fNonNan, fPriceVolume)
499
500 for sAdd in lsValid:
501 if sAdd in lsCheckStocks:
502 lsRetStocks.append( sAdd )
503
504 return sorted(lsRetStocks)
505