Textual Analysis on S&P500 Companies

Bag of words and Doc2Vec on company business descriptions

Data
S&P 500 membership and identifiers from CRSP (msp500list_v2, msf_v2, stksecurityinfohist) CIZ format; business descriptions from Compustat and Capital IQ
Sample
S&P 500 companies as of December 2020
Methods
Text cleaning with NLTK and spaCy, bag of words with scikit-learn, Doc2Vec with gensim
Author
Qingyi (Freda) Song Drechsler

Download notebook (.ipynb)

This set of code demonstrates how to conduct textual analysis using S&P500 companies’ business description. It covers some basic textual analysis steps, including:

part 0: import packages

This code uses many NLP related packages, such as NLTK for stop words, Gensim and SpaCy for various preparation, vectorization tasks, and scikit-learn for vectorization and similarity calculation. Of course, there are many other ways of handling this task, and this code just serves as a toy model to demonstrate the concept.

#########################################
# Textual Analysis on S&P 500 Companies #
# CRSP CIZ Format                       #
# Qingyi (Freda) Song Drechsler         #
# Date: April 2021                      #
# Updated: October 2026                 #
#########################################

# Import packages

# wrds
import wrds

# common packages
import os
import re
import pandas as pd
import numpy as np
import pickle as pkl
from pprint import pprint
from collections import OrderedDict
from collections import defaultdict
import multiprocessing 

import nltk
from nltk.corpus import stopwords
stop_words = stopwords.words("english")

# gensim
import gensim
import gensim.corpora as corpora
from gensim.utils import simple_preprocess
from gensim.models import CoherenceModel
from gensim.test.utils import common_corpus, common_dictionary
from gensim.models.doc2vec import Doc2Vec, TaggedDocument

# spacy
import spacy
nlp = spacy.load('en_core_web_sm', disable=['parser', 'ner'])

# sklearn
from sklearn.metrics.pairwise import cosine_similarity
from sklearn.metrics import accuracy_score
from sklearn.feature_extraction.text import TfidfVectorizer
from sklearn.feature_extraction.text import CountVectorizer
from sklearn.decomposition import LatentDirichletAllocation

# Visualization
import matplotlib.pyplot as plt
%matplotlib inline

part 1: build S&P500 index constituents

Let’s first establish connection with WRDS server and extract the constituents information of S&P500 Index.

###################
# Connect to WRDS #
###################
conn=wrds.Connection()


### Get S&P500 Index Membership from CRSP
### I opt for the monthly frequency of the data,
### but one can choose to work with crsp.dsp500list_v2
### if more precise date range is needed.

sp500 = conn.raw_sql("""
    select distinct a.permno, a.mbrstartdt, a.mbrenddt,
           b.mthcaldt, b.mthret
    from crsp_m_indexes.msp500list_v2 as a,
         crsp_m_stock.msf_v2 as b
    where a.permno=b.permno
    and b.mthcaldt >= a.mbrstartdt and b.mthcaldt <= a.mbrenddt
    and b.mthcaldt between '01/01/2000' and '06/30/2026'
    order by b.mthcaldt, a.permno
    """, date_cols=['mbrstartdt', 'mbrenddt', 'mthcaldt'])
Loading library list...
Done

Prior to 2020, Compustat used to provide index constituents data through WRDS, however, this data component was discontinued by Compustat as of October 2020, and researchers can now access S&P500 Index constituents through CRSP. However, many of the other index constituents information are no longer available due to coverage difference between CRSP and Compustat. If you as a researcher find this change inconvenient for your research projects, please reach out to S&P support support.datafeed.mi@spglobal.com directly.

Next, let’s add various other company identifiers (e.g. ticker) by linking the sp500 dataframe with crsp.msenames.

### Add Other Company Identifiers from CRSP.STKSECURITYINFOHIST
### - You don't need this step if only PERMNO is required
### - This step aims to add TICKER, CUSIP, company name,
###   share type, primary exchange, SIC code and etc.

stkinfo = conn.raw_sql("""
    select distinct permno, issuernm, cusip, ticker,
           sharetype, securitytype, securitysubtype,
           primaryexch, siccd, secinfostartdt, secinfoenddt
    from crsp_m_stock.stksecurityinfohist
    """, date_cols=['secinfostartdt', 'secinfoenddt'])

# Merge with SP500 data
sp500_full = pd.merge(sp500, stkinfo, how='left', on='permno')

# Impose the date range restrictions
sp500_full = sp500_full.loc[
    (sp500_full.mthcaldt >= sp500_full.secinfostartdt)
    & (sp500_full.mthcaldt <= sp500_full.secinfoenddt)]

As the corpus of this code comes from CIQ’s company business description, we need to add company identifier used in the CIQ eco system, GVKEY.

### Add Compustat Identifiers
### - Link with Compustat's GVKEY and IID if need to work with
###   fundamental data
### - Linkage is done through crsp.ccmxpf_linktable

ccm=conn.raw_sql("""
    select distinct gvkey, liid as iid, lpermno as permno,
           linktype, linkprim, linkdt, linkenddt
    from crsp_m_ccm.ccmxpf_linktable
    where substr(linktype,1,1)='L'
    and (linkprim ='C' or linkprim='P')
    """, date_cols=['linkdt', 'linkenddt'])

# if linkenddt is missing then set to today date
ccm['linkenddt']=ccm['linkenddt'].fillna(pd.to_datetime('today'))

# Merge the CCM data with S&P500 data
# First just link by matching PERMNO
sp500ccm = pd.merge(sp500_full, ccm, how='left', on=['permno'])

# Then set link date bounds
sp500ccm = sp500ccm.loc[(sp500ccm['mthcaldt']>=sp500ccm['linkdt'])\
                        &(sp500ccm['mthcaldt']<=sp500ccm['linkenddt'])]

# Rearrange columns for final output

sp500ccm = sp500ccm[['mthcaldt', 'permno', 'issuernm', 'cusip',
                     'ticker', 'sharetype', 'primaryexch', 'siccd',
                     'gvkey', 'iid', 'mbrstartdt', 'mbrenddt',
                     'mthret']]

Optional step of linking with CIK in case you want to connect with SEC filings for your analysis:

### Add CIKs and Link with SEC Index Files using CIK

names = conn.raw_sql("""
    select gvkey, cik, sic, naics, gind, gsubind
    from comp.names
    """)

# Merge sp500 constituents table with names table
sp500cik = pd.merge(sp500ccm, names, on='gvkey',  how='left')
sp500cik.head()
mthcaldt permno issuernm cusip ticker sharetype primaryexch siccd gvkey iid mbrstartdt mbrenddt mthret cik sic naics gind gsubind
0 2000-01-31 10078 SUN MICROSYSTEMS INC 86681010 SUNW NS Q 3570 012136 01 1992-08-20 2010-01-28 0.014528 0000709519 3571 334111 452020 45202010
1 2000-01-31 10104 ORACLE CORP 68389X10 ORCL NS Q 7370 012142 01 1989-08-03 2026-08-31 -0.108477 0001341439 7372 513210 451030 45103020
2 2000-01-31 10107 MICROSOFT CORP 59491810 MSFT NS Q 7370 012141 01 1994-06-07 2026-08-31 -0.161670 0000789019 7372 513210 451030 45103020
3 2000-01-31 10138 T ROWE PRICE ASSOC INC 74147710 TROW NS Q 6210 012138 01 1999-10-13 2026-08-31 0.052454 0001113169 6282 523940 402030 40203010
4 2000-01-31 10145 HONEYWELL INTERNATIONAL INC 43851610 HON NS N 3714 001300 01 1925-12-31 2026-08-31 -0.167931 0000773840 9997 999977 201050 20105010

part 2: read in business description

S&P offers various databases under its umbrella, including most commonly used Compustat (COMP) and Capital IQ (CIQ). This block of code reads in the long format of business description from CIQ and short version from COMP. We will be using the long version as corpus for this project.

### Extract business description

# short description - "busdesc" from comp.company 
compbus = conn.raw_sql("""
    select gvkey, conm, cik, busdesc from comp.company
    """)

# long description - "businessdescription" from CIQ
# link with GVKEY using CompanyID
ciqbus = conn.raw_sql("""
    select a.*, b.gvkey, b.startdate, b.enddate
    from ciq.ciqbusinessdescription as a,
    ciq.wrds_gvkey as b
    where a.companyid = b.companyid
    """, date_cols = ['startdate', 'enddate'])

Occasionally, a GVKEY can be mapped with multiple CompanyIDs (identifiers used in CIQ system). Here is one example:

# GVKEY 003413 maps to more than one CIQ CompanyID
ciqbus.loc[ciqbus.gvkey == '003413']
companyid businessdescription gvkey startdate enddate
6755 263295.0 Consolidated Edison, Inc., through its subsidi... 003413 1998-01-02 NaT
43234 6020983.0 Consolidated Edison Company of New York, Inc. ... 003413 NaT 1998-01-01

To avoid double counting, for a given GVKEY, if more than one companyid is mapped, I pick the larger companyid as it is probably a more recent link. This will yield a one-to-one mapping between GVKEY and CompanyID.

# Choose record with max CompanyID
maxid = ciqbus.groupby('gvkey')['companyid'].max().reset_index()

ciqbus_nodup = pd.merge(ciqbus, maxid, how = 'inner',
                        on = ['gvkey', 'companyid'])
ciqbus_nodup.loc[ciqbus_nodup.gvkey == '003413']
companyid businessdescription gvkey startdate enddate
39989 6020983.0 Consolidated Edison Company of New York, Inc. ... 003413 NaT 1998-01-01

Let’s inspect the short and long version of business description of the same company, using AAI Corporation as an example.

### Merge business description from COMP and CIQ
# join by GVKEY
busdesc = pd.merge(compbus, ciqbus_nodup, how = 'inner', on = 'gvkey')

# Print one short and long version of the business description
# for inspection

# AAI Corporation (GVKEY 001002); selected by GVKEY,
# as row order can change
aai = busdesc.loc[busdesc.gvkey == '001002'].iloc[0]
print('\nShort business description from COMP:\n', aai.busdesc )
print('\nLong business description from CIQ:\n',  aai.businessdescription )

Short business description from COMP:
 Textron Systems Corporation designs, develops, manufactures, and supports aerospace and defense technologies. It offers unmanned aircraft systems, micro air vehicles, remote video terminals, ground control stations, and portable ground control stations.

Long business description from CIQ:
 Textron Systems Corporation designs, develops, manufactures, and supports aerospace and defense technologies. It offers unmanned aircraft systems, micro air vehicles, remote video terminals, ground control stations, and portable ground control stations. The company also designs and develops training solutions, such as sonar, satellite, and aircraft trainers. In addition, it offers simulation-based maintenance training devices, modular simulated aircraft maintenance trainers, sentry airborne warning and control systems, joint surveillance target attack radar systems, electronic warfare trainers, embedded on-board trainers, and air defense trainers. Further, the company provides flight line, depot, and product line electronic test equipment; and acoustic detection systems, aviation ground support equipment, and armament and mechanical systems. Furthermore, it offers operations, maintenance, logistics, and supply chain services; and engineering and modification, training and simulation, and depot support services. The company offers its products for intelligence gathering, surveillance, target acquisition, reconnaissance, observation, and battle damage assessment applications. Textron Systems Corporation has strategic partnership with Kaman. Textron Systems Corporation was founded in 1950 and is based in Hunt Valley, Maryland. Textron Systems Corporation operates as a subsidiary of Textron Inc.

part 3: clean and prepare corpus

This is a standard cleaning process that can be applied to other text bodies, user can modify the procedure to suite specific textual cleaning need. It does the following common cleaning tasks:

  • convert text to lower case
  • keep only alphabetic characters (remove numbers and special characters)
  • remove stop words
  • form bigrams
  • lemmatize words

Let’s use companies in the S&P500 Index as of 2020 as our universe:

# Use SP500 Companies in 2020 as universe
univ = sp500cik.loc[sp500cik.mthcaldt=='2020-12-31'][
    ['mthcaldt', 'permno', 'issuernm', 'cusip', 'gvkey', 'iid',
     'cik', 'ticker', 'sic', 'naics']]

# Merge to get business description
sp500busdesc = pd.merge(univ, busdesc[['gvkey', 'companyid', 'busdesc',
                                       'businessdescription']],
                        how = 'inner', on = ['gvkey'])
sp500busdesc.head()
mthcaldt permno issuernm cusip gvkey iid cik ticker sic naics companyid busdesc businessdescription
0 2020-12-31 10104 ORACLE CORP 68389X10 012142 01 0001341439 ORCL 7372 513210 22247.0 Oracle Corporation offers products and service... Oracle Corporation offers products and service...
1 2020-12-31 10107 MICROSOFT CORP 59491810 012141 01 0000789019 MSFT 7372 513210 21835.0 Microsoft Corporation, a technology company, d... Microsoft Corporation, a technology company, d...
2 2020-12-31 10138 T ROWE PRICE GROUP INC 74144T10 012138 01 0001113169 TROW 6282 523940 23217.0 T. Rowe Price Group, Inc. is a publicly owned ... T. Rowe Price Group, Inc. is a publicly owned ...
3 2020-12-31 10145 HONEYWELL INTERNATIONAL INC 43851610 001300 01 0000773840 HON 9997 999977 1340740.0 Honeywell International Inc. engages in the in... Honeywell International Inc. engages in the in...
4 2020-12-31 10516 ARCHER DANIELS MIDLAND CO 03948310 001722 01 0000007084 ADM 2070 31122 251704.0 Archer-Daniels-Midland Company provides human ... Archer-Daniels-Midland Company provides human ...

Start Cleaning. First convert the busdesc (short) and businessdescription (long) column into lists.

# Convert both short and long busindess description
# (busdesc, businessdescription) to list
busdesc_lst = sp500busdesc.busdesc.tolist()
busdesclong_lst = sp500busdesc.businessdescription.tolist()

# For this exercise, I use long business description as corpus
text_feed = busdesclong_lst

len(busdesc_lst), len(busdesclong_lst)
(501, 501)

Below are several functions for textual cleaning tasks:

  • text_clean() - convert all texts to lower case and remove non-letter components
# convert all text to lower case and remove non alphabetic components
def text_clean(raw_text):
    # 1. convert to lower case
    txt1 = raw_text.lower().strip()
    
    # 2. remove non-letters 
    txt2 = re.sub("[^a-zA-Z]", " ", txt1).split(' ')
    return (txt2)
  • sent_to_words() - convert sentences to list
# convert sentences to list
def sent_to_words(sentences):
    for sent in sentences:
        yield(text_clean(sent))
  • remove_stopwords()- take out the stop words defined by NLTK
# Remove Stop Words
def remove_stopwords(texts):
    nostop = [[word for word in doc if word not in stop_words]
              for doc in texts]
    return nostop
  • make_bigrams()- forms bigram if needed
# form bigrams
def make_bigrams(texts):
    #bigram = gensim.models.Phrases(texts, min_count=1, threshold=10)
    bigram = gensim.models.Phrases(texts)
    bigram_mod = gensim.models.phrases.Phraser(bigram)
    return [bigram_mod[doc] for doc in texts]
  • lemmatization_str()- lemmatize the words and convert the output to list of strings
# lemmatization with output of list of strings to fit into countvectorizer
def lemmatization_str(texts):
    texts_out = []
    for sent in texts:
        doc = nlp(" ".join(sent))
        texts_out.append(' '.join([token.lemma_ for token in doc]))
    return texts_out

Now let’s apply these functions to the input list:

### Apply the raw text into the cleaning process
# Input text: text_feed
# Output text: text_prep

print('Start cleaning ... ')
print('Remove non-alphabetic characters and convert to lower case ... ')
text_nonum = list(sent_to_words(text_feed))

print('Remove stop words ... ')
text_nostops = remove_stopwords(text_nonum)

print('Form bigrams ... ')
text_bigrams = make_bigrams(text_nostops)

print('Lemmatize words ... ')
text_lemmatized = lemmatization_str(text_bigrams)

print('Cleaning done. text_prep is the cleaned output.')
text_prep = text_lemmatized
len(text_prep)
Start cleaning ... 
Remove non-alphabetic characters and convert to lower case ... 
Remove stop words ... 
Form bigrams ... 
Lemmatize words ... 
Cleaning done. text_prep is the cleaned output.
501

The cleaning process above should convert the input text to all lower case, remove non letters components, form bigrams and lemmatize the words. Below is a comparison of input texts versus the cleaned text, using Ralph Lauren’s business description as an example:

Raw input text:

# Raw input text: Ralph Lauren's long business description
rl_pos = sp500busdesc.loc[sp500busdesc.ticker == 'RL'].index.values[0]
print(text_feed[rl_pos])
Ralph Lauren Corporation designs, markets, and distributes lifestyle products in North America, Europe, Asia, and internationally. The company offers apparel, including a range of men’s, women’s, and children’s clothing; footwear and accessories, which comprise casual shoes, dress shoes, boots, sneakers, sandals, eyewear, watches, fashion and fine jewelry, scarves, hats, gloves, and umbrellas, as well as leather goods comprising handbags, luggage, small leather goods, and belts; home products, such as bed and bath lines, furniture, fabric and wall coverings, lighting, dining, floor coverings, decorative accessories, and giftware; and fragrances. The company sells apparel and accessories under the Ralph Lauren Collection, Ralph Lauren Purple Label, Polo Ralph Lauren, Double RL, Lauren Ralph Lauren, RLX Ralph Lauren, Polo Ralph Lauren Children, and Chaps brands; women's fragrances under the Ralph Lauren Collection, Woman by Ralph Lauren, Romance Collection, and Ralph Collection brand names; and men's fragrances under the Ralph's Club, Purple Label, Polo Blue, Polo Red, Polo Green, Polo Black, Polo 67, Safari, Polo Sport, and Big Pony Men's brand names. Its restaurant collection includes The Polo Bar in New York City; RL Restaurant in Chicago; Ralph's in Paris; The Bar at Ralph Lauren located in Milan; Ralph's Bar located in Chengdu, China; and Ralph's Coffee concept. The company sells its products to department stores, specialty stores, and golf and pro shops, as well as directly to consumers through its retail stores, concession-based shop-within-shops, and its digital commerce sites. It operates retail stores and concession-based shop-within-shops; and operates Ralph Lauren stores and shops through licensing partners. Ralph Lauren Corporation was founded in 1967 and is based in New York, New York.

Cleaned output text:

# Cleaned output text
print(text_prep[rl_pos])
ralph_lauren corporation design market distribute lifestyle product north_america europe asia internationally company_offers apparel include range man woman child clothing footwear accessory comprise casual shoe dress shoe boot sneaker sandal eyewear watch fashion fine jewelry scarf hat glove umbrella well leather good comprise handbag luggage small leather good belt home product bed bath line furniture fabric wall covering lighting dining floor covering decorative accessory giftware fragrance company sell apparel accessory ralph_lauren collection ralph_lauren purple label polo ralph_lauren double rl lauren ralph_lauren rlx ralph_lauren polo ralph_lauren child chap brand woman fragrance ralph_lauren collection woman ralph_lauren romance collection ralph collection brand_name man fragrance ralph club purple label polo blue polo red polo green polo black polo safari polo sport big pony man brand_names restaurant collection include polo bar new_york city rl restaurant chicago ralph paris bar ralph_lauren locate milan ralph bar locate chengdu china ralph coffee concept company sells_product department_stores specialty store golf pro shop well directly_consumer retail_store concession base shop within shop digital commerce_sites operate retail_store concession base shop within shop operate ralph_lauren store shop license partner ralph_lauren corporation_founde base new_york new_york

part 4: form bag of words

This block of code uses the vectorization tools offered in scikit-learn, and forms bag-of-words (BOW) for each corpus.

# Vectorize the vocab
vect = CountVectorizer()

# Form vector for corpus
vect_array = vect.fit_transform(text_prep).toarray()

### Create BOW Dictionary

intext = text_prep

# Create an empty dictionary to store the similarity score
bow_dict = {}

for row_id, text in enumerate(intext):
    # loop through all records
    frequency = defaultdict(int)
    
    for token in text.split(' '):
        #store tokens and count in frequency dict
        frequency[token] += 1
    
    # store frequency dict to bow_dict
    bow_dict[row_id] = frequency

# Store BOW into df
bow_df = pd.DataFrame()
for item, lst in list(bow_dict.items()):
    _df = pd.DataFrame(lst.keys(), lst.values())
    _df['comp_id'] = item
    bow_df = pd.concat([bow_df, _df])

bow_df = bow_df.reset_index().rename(columns={'index': 'freq',
                                              0: 'word'})

# slice out SP500 company ids for linking
sp500id = sp500busdesc[['permno', 'gvkey', 'iid', 'ticker']]\
    .reset_index().rename(columns={'index': 'comp_id'})

# Add PERMNO GVKEY IID info
bow_df = pd.merge(bow_df, sp500id, how = 'left', on = ['comp_id'])

The BOW dataframe contains the following information:

  • word: focal word
  • freq: number of times a word appears in that document
  • company identifiers

Below is a subset of BOW for the company Apple. For instance, the word “apple” appears 15 times, and the word “service” appears 6 times in Apple’s business description obtained from CIQ.

# Most frequent words in Apple's business description
bow_df.loc[bow_df.ticker == 'AAPL'][['word', 'freq', 'permno',
                                     'gvkey', 'iid', 'ticker']]\
    .sort_values(by='freq', ascending=False).head(10)
word freq permno gvkey iid ticker
5445 apple 15 14593 001690 01 AAPL
5481 service 6 14593 001690 01 AAPL
5476 third_party 4 14593 001690 01 AAPL
5446 inc 4 14593 001690 01 AAPL
5455 accessory 3 14593 001690 01 AAPL
5487 store 3 14593 001690 01 AAPL
5459 line 3 14593 001690 01 AAPL
5474 well 3 14593 001690 01 AAPL
5452 computer 3 14593 001690 01 AAPL
5504 subscription 3 14593 001690 01 AAPL

And out of curiosity, I summarize below the top 10 most frequently appeared words in all S&P500 companies’ business description. Researchers might want to customize the generic stop words defined by NLTK to include certain content specific stop words, such as “company” or “inc” in this case.

### Examine word distribution across all company business
### descriptions (documents)

totalcnt = bow_df.groupby('word')['freq'].sum().reset_index()\
    .sort_values(by=['freq'], ascending = False)

# Report top 10 most frequent words
totalcnt.head(10)
word freq
7248 service 1262
6354 product 1174
1689 company 932
7489 solution 705
4007 inc 666
7894 system 488
4012 include 484
8794 well 445
7190 segment 436
6420 provide 411

part 5: similarity based on Bow

This section relies on the vectors created based on BOW and calculates the pairwise similarity score among companies in this universe. Before diving into the full sample, let’s inspect this approach using individual companies: Disney, Microsoft, Wells Fargo and Citi.

### Illustrating using several companies: 

# Find position of the document in the list for each company
dis_pos  = sp500busdesc.loc[sp500busdesc.ticker =='DIS'].index.values[0]
msft_pos = sp500busdesc.loc[sp500busdesc.ticker =='MSFT'].index.values[0]
wfc_pos  = sp500busdesc.loc[sp500busdesc.ticker =='WFC'].index.values[0]
c_pos    = sp500busdesc.loc[sp500busdesc.ticker =='C'].index.values[0]

### Report pairwise similarity score among the select firms
# use cosine_similarity function
print('Pairwise similarity score:')
def bow_pair(i, j):
    return cosine_similarity([vect_array[i]], [vect_array[j]])[0, 0]

print('Microsoft vs Disney:',
      "{:.2f}".format(bow_pair(dis_pos, msft_pos)))
print('Disney vs Wells Fargo:',
      "{:.2f}".format(bow_pair(dis_pos, wfc_pos)))
print('Wells Fargo vs Citi:',   "{:.2f}".format(bow_pair(c_pos, wfc_pos)))
Pairwise similarity score:
Microsoft vs Disney: 0.15
Disney vs Wells Fargo: 0.14
Wells Fargo vs Citi: 0.61

And the pairwise similarity score among these companies are reported below. As anticipated, the two financial institutions, Wells Fargo and Citi, display relatively high similarity score based on BOW formed from their respective business description.

Now let’s apply the similarity score calculation on the entire sample:

### Similarity score across all companies

# list containing all tickers
ticker_lst = sp500busdesc.ticker.tolist()

# empty list to store the bow similarity scores
bow_score = []

# input vector
invect = vect_array

# function to loop through company pair (i,j)
def bow_simscore(inlist):
    for i in range(len(invect)):
        for j in range(len(invect)):
            
            _tmp = {'comp1': ticker_lst[i],  # ticker for company i
                    'comp2': ticker_lst[j],  # ticker for company j
                    # similarity score for pair (i,j)
                    'bow_score': cosine_similarity([invect[i]],
                                                   [invect[j]])[0, 0]}
            bow_score.append(_tmp)
    return bow_score    

### Apply over all companies
print ('Calculating pairwise similarity score among all firms ... ')
bow_score = bow_simscore(invect)
print ('Done.')
Calculating pairwise similarity score among all firms ... 
Done.

Tidy up the output and sort the results in descending order to find the most similar company for each focal company.

### Tidying up the output

# Convert to dataframe 
bow_score_df = pd.DataFrame(bow_score)

# take out comparison against itself (i,i), which is always 1 by nature
# then drop the (i,j) vs (j, i) duplicates
bow_score_df = bow_score_df.loc[bow_score_df.comp1 != bow_score_df.comp2]\
    .drop_duplicates()

# Sort by comp1 then descending based on bow_score
# to find most similar companies
bow_score_df = bow_score_df.sort_values(by=['comp1', 'bow_score'],
                                        ascending = [True, False])
bow_score_df.head()
comp1 comp2 bow_score
210744 A ROP 0.420425
210683 A TMO 0.411720
210792 A WAT 0.361186
210488 A KEYS 0.360261
210638 A DHR 0.356645

Then for each of the S&P 500 companies, find its most similar peer based on BOW similarity score:

# For each company, keep the company with the highest cosine score
similar_comp_bow = bow_score_df.loc[
    bow_score_df.groupby('comp1')['bow_score'].idxmax()]

# merge back company identifier info
_ids = sp500busdesc[['issuernm', 'ticker', 'sic', 'naics']]
similar_comp_bow = pd.merge(similar_comp_bow, _ids, how = 'left',
                            left_on = 'comp1', right_on = 'ticker')
similar_comp_bow = pd.merge(similar_comp_bow, _ids, how = 'left',
                            left_on = 'comp2', right_on = 'ticker')
similar_comp_bow = similar_comp_bow.rename(columns={
    'comp2': 'comp_bow',
    'issuernm_x':'issuernm1', 'issuernm_y':'issuernm_bow',
    'sic_x': 'sic1', 'sic_y':'sic_bow',
    'naics_x':'naics1', 'naics_y':'naics_bow'})\
    .drop(columns=['ticker_x', 'ticker_y'])
similar_comp_bow = similar_comp_bow[
    ['comp1', 'issuernm1', 'sic1', 'naics1', 'comp_bow',
     'issuernm_bow', 'sic_bow', 'naics_bow', 'bow_score']]\
    .drop_duplicates()

I report below the most similar company for some of the familiar tickers. Based on BOW pairwise similarity score, some company has clear peers (e.g. GM vs F and or V vs MA), others report most similar peer somehow different from their dominant industry (e.g. FB vs TWTR or AMZN vs ADBE). The latter case often happens among companies that expand in multiple fields of operation.

# Lookup most similar company by ticker
company_ticker = ['NCLH', 'MSFT', 'WFC', 'XOM', 
                  'AMZN', 'C', 'DIS', 'V', 'JNJ',
                 'WMT', 'UAL', 'FB', 'CBOE', 
                 'GM', 'ABT', 'MMM', 'BMY']

sample_bow = similar_comp_bow.loc[
    similar_comp_bow.comp1.isin(company_ticker)]
sample_bow.style.format({'bow_score': '{:.2f}'})
  comp1 issuernm1 sic1 naics1 comp_bow issuernm_bow sic_bow naics_bow bow_score
7 ABT ABBOTT LABORATORIES 3845 334510 BDX BECTON DICKINSON & CO 3841 339112 0.45
35 AMZN AMAZON COM INC 5961 455219 ADBE ADOBE INC 7370 518210 0.39
68 BMY BRISTOL MYERS SQUIBB CO 2834 325412 AMGN AMGEN INC 2836 325414 0.48
74 C CITIGROUP INC 6199 522110 USB U S BANCORP DEL 6020 522110 0.64
80 CBOE C B O E GLOBAL MARKETS INC 6200 523210 CME C M E GROUP INC 6200 523210 0.42
133 DIS DISNEY WALT CO 4888 516120 FOXA FOX CORP 4833 516120 0.38
175 FB FACEBOOK INC 7370 519290 TWTR TWITTER INC 7370 519130 0.25
200 GM GENERAL MOTORS CO 3711 336110 F FORD MOTOR CO DEL 3711 336110 0.52
253 JNJ JOHNSON & JOHNSON 2834 325412 ZBH ZIMMER BIOMET HOLDINGS INC 3842 339113 0.41
306 MMM 3M CO 9997 999977 NWL NEWELL BRANDS INC 3630 3352 0.46
315 MSFT MICROSOFT CORP 7372 513210 FFIV F 5 NETWORKS INC 7373 541512 0.45
321 NCLH NORWEGIAN CRUISE LINE HLDGS LTD 4400 483112 CCL CARNIVAL CORP 4400 483112 0.54
449 UAL UNITED AIRLINES HOLDINGS INC 4512 481111 EXPE EXPEDIA GROUP INC 4700 561510 0.37
459 V VISA INC 6099 522320 MA MASTERCARD INC 7389 522320 0.45
480 WFC WELLS FARGO & CO NEW 6020 522110 USB U S BANCORP DEL 6020 522110 0.66
485 WMT WALMART INC 5331 455219 TGT TARGET CORP 5331 455219 0.35
494 XOM EXXON MOBIL CORP 2911 324110 PSX PHILLIPS 66 2911 324110 0.52

part 6: similarity based on doc2vec

This is an alternative approach to calculate similarity score, and it relies on Doc2Vec as part of the Gensim package. Please note that the sample size of this toy project is probably too small but I nevertheless include it here for demonstration purpose.

# Uses text_nonum as input as doc2vec approach doesn't require
# cleaning of stop words or lemmatization

# Tagging documents
tagged_text = [TaggedDocument(words=_d, tags=[str(i)])
               for i, _d in enumerate(text_nonum)]

# Now train the data
#cores = multiprocessing.cpu_count()
assert gensim.models.doc2vec.FAST_VERSION > -1 # accelerate the process

# Build model
# fixed seed and a single worker make the results reproducible
model = Doc2Vec(dm=1, vector_size=40, min_count = 1, seed=2021, workers=1)
model.build_vocab(tagged_text)

for epoch in range(10):
    model.train(tagged_text, epochs=model.epochs,
                total_examples=model.corpus_count)
    print("Epoch #{} is trained.".format(epoch+1))
Epoch #1 is trained.
Epoch #2 is trained.
Epoch #3 is trained.
Epoch #4 is trained.
Epoch #5 is trained.
Epoch #6 is trained.
Epoch #7 is trained.
Epoch #8 is trained.
Epoch #9 is trained.
Epoch #10 is trained.

Again let’s first inspect the pairwise similarity score among the four individual companies we examined before.

### Report pairwise similarity score among the select firms
# use cosine_similarity function
print('Pairwise similarity score:')
def d2v_pair(i, j):
    return model.wv.n_similarity(text_nostops[i], text_nostops[j])

print('Microsoft vs Disney:',
      "{:.2f}".format(d2v_pair(msft_pos, dis_pos)))
print('Disney vs Wells Fargo:',
      "{:.2f}".format(d2v_pair(dis_pos, wfc_pos)))
print('Wells Fargo vs Citi:',   "{:.2f}".format(d2v_pair(wfc_pos, c_pos)))
Pairwise similarity score:
Microsoft vs Disney: 0.79
Disney vs Wells Fargo: 0.63
Wells Fargo vs Citi: 0.96

Let’s now expand the analysis to the entire sample of all 500 companies:

# Calculate d2v based similarity score for the entire sample
d2v_score = []

def d2v_simscore(inlist):
    for i in range(len(inlist)):
        for j in range(len(inlist)):
            
            _tmp = {'comp1': ticker_lst[i], 
                    'comp2': ticker_lst[j],
                    'd2v_score':  model.wv.n_similarity(inlist[i],
                                                         inlist[j]) }
            d2v_score.append(_tmp)
    return d2v_score    

# apply over all companies
d2v_score = d2v_simscore(text_nostops)

# Convert to dataframe 
d2v_score_df = pd.DataFrame(d2v_score)

# take out comparison against itself (i,i)
# then drop the (i,j) vs (j, i) duplicates
d2v_score_df = d2v_score_df.loc[d2v_score_df.comp1 != d2v_score_df.comp2]\
    .drop_duplicates()

# Sort by comp1 then descending based on bow_score
d2v_score_df = d2v_score_df.sort_values(by=['comp1', 'd2v_score'],
                                        ascending = [True, False])

And then pick out the most similar company:

# For each company, keep the most similar company
similar_comp_d2v = d2v_score_df.loc[
    d2v_score_df.groupby('comp1')['d2v_score'].idxmax()]

similar_comp_d2v = pd.merge(similar_comp_d2v, _ids, how = 'left',
                            left_on = 'comp1', right_on = 'ticker')
similar_comp_d2v = pd.merge(similar_comp_d2v, _ids, how = 'left',
                            left_on = 'comp2', right_on = 'ticker')
similar_comp_d2v = similar_comp_d2v.rename(columns={
    'comp2': 'comp_d2v',
    'issuernm_x':'issuernm1', 'issuernm_y':'issuernm_d2v',
    'sic_x': 'sic1', 'sic_y':'sic_d2v',
    'naics_x':'naics1', 'naics_y':'naics_d2v'})\
    .drop(columns=['ticker_x', 'ticker_y'])
similar_comp_d2v = similar_comp_d2v[
    ['comp1', 'issuernm1', 'sic1', 'naics1', 'comp_d2v',
     'issuernm_d2v', 'sic_d2v', 'naics_d2v', 'd2v_score']]\
    .drop_duplicates()

Again, let’s inspect the same group of companies before:

sample_d2v = similar_comp_d2v.loc[
    similar_comp_d2v.comp1.isin(company_ticker)]
sample_d2v.style.format({'d2v_score': '{:.2f}'})
  comp1 issuernm1 sic1 naics1 comp_d2v issuernm_d2v sic_d2v naics_d2v d2v_score
7 ABT ABBOTT LABORATORIES 3845 334510 BDX BECTON DICKINSON & CO 3841 339112 0.95
35 AMZN AMAZON COM INC 5961 455219 ADBE ADOBE INC 7370 518210 0.94
68 BMY BRISTOL MYERS SQUIBB CO 2834 325412 ABBV ABBVIE INC 2836 325414 0.96
74 C CITIGROUP INC 6199 522110 JPM JPMORGAN CHASE & CO 6020 522110 0.97
80 CBOE C B O E GLOBAL MARKETS INC 6200 523210 ICE INTERCONTINENTALEXCHANGE GRP INC 6200 523999 0.93
132 DIS DISNEY WALT CO 4888 516120 VIAC VIACOMCBS INC 4833 516120 0.96
174 FB FACEBOOK INC 7370 519290 EA ELECTRONIC ARTS INC 7372 513210 0.91
199 GM GENERAL MOTORS CO 3711 336110 F FORD MOTOR CO DEL 3711 336110 0.97
252 JNJ JOHNSON & JOHNSON 2834 325412 SYK STRYKER CORP 3841 339112 0.94
306 MMM 3M CO 9997 999977 ITW ILLINOIS TOOL WORKS INC 3560 33399 0.96
315 MSFT MICROSOFT CORP 7372 513210 CSCO CISCO SYSTEMS INC 3576 334118 0.94
321 NCLH NORWEGIAN CRUISE LINE HLDGS LTD 4400 483112 CCL CARNIVAL CORP 4400 483112 0.94
449 UAL UNITED AIRLINES HOLDINGS INC 4512 481111 BKNG BOOKING HOLDINGS INC 4700 561510 0.94
459 V VISA INC 6099 522320 MA MASTERCARD INC 7389 522320 0.96
480 WFC WELLS FARGO & CO NEW 6020 522110 KEY KEYCORP NEW 6020 522110 0.98
485 WMT WALMART INC 5331 455219 NKE NIKE INC 3021 316210 0.94
494 XOM EXXON MOBIL CORP 2911 324110 CVX CHEVRON CORP NEW 2911 324110 0.96

And finally, let’s put the BOW and Doc2Vec results side by side for comparison. For some instances, both approaches yield the same result for the focal company (e.g. ABT vs BDX, GM vs F, NCLH vs CCL, V vs MA), others return different candidates but generally belong to the same industry.

# BOW and Doc2Vec results side by side
compare = pd.merge(sample_bow[['comp1', 'issuernm1',
                               'comp_bow', 'bow_score']],
                   sample_d2v[['comp1', 'comp_d2v', 'd2v_score']],
                   how='inner', on='comp1')
compare.style.format({'bow_score': '{:.2f}', 'd2v_score': '{:.2f}'})
  comp1 issuernm1 comp_bow bow_score comp_d2v d2v_score
0 ABT ABBOTT LABORATORIES BDX 0.45 BDX 0.95
1 AMZN AMAZON COM INC ADBE 0.39 ADBE 0.94
2 BMY BRISTOL MYERS SQUIBB CO AMGN 0.48 ABBV 0.96
3 C CITIGROUP INC USB 0.64 JPM 0.97
4 CBOE C B O E GLOBAL MARKETS INC CME 0.42 ICE 0.93
5 DIS DISNEY WALT CO FOXA 0.38 VIAC 0.96
6 FB FACEBOOK INC TWTR 0.25 EA 0.91
7 GM GENERAL MOTORS CO F 0.52 F 0.97
8 JNJ JOHNSON & JOHNSON ZBH 0.41 SYK 0.94
9 MMM 3M CO NWL 0.46 ITW 0.96
10 MSFT MICROSOFT CORP FFIV 0.45 CSCO 0.94
11 NCLH NORWEGIAN CRUISE LINE HLDGS LTD CCL 0.54 CCL 0.94
12 UAL UNITED AIRLINES HOLDINGS INC EXPE 0.37 BKNG 0.94
13 V VISA INC MA 0.45 MA 0.96
14 WFC WELLS FARGO & CO NEW USB 0.66 KEY 0.98
15 WMT WALMART INC TGT 0.35 NKE 0.94
16 XOM EXXON MOBIL CORP PSX 0.52 CVX 0.96