{mapBliss} is a package initially developed in my attempt to copy the offerings made available by a site (since defunct) called Atlas.co. Recently I was talking about the package in a group chat and someone reached out to me sharing that he traveled to 43 states in the United States and shared with me a list of his journey.

These are the visuals that were made:

In this blog, I share a new dataset available in {mapBliss} – us_travels. This dataset provides an interesting example for plotting more lengthy and complex adventures with {mapBliss}.

Loading the us_travels dataset

To load and preview the us_travels dataset, run the following code:

# Install the mapBliss package from GitHub if you haven't already
# devtools::install_packages("benyamindsmith/mapBliss")

# load the mapBliss library
library(mapBliss)
# load the us_travels dataset
data("us_travels")

# Preview the data
us_travels

# A tibble: 85 × 2
   location                       how  
   <chr>                          <chr>
 1 Mobile, Alabama, USA           car  
 2 Montgomery, Alabama, USA       car  
 3 Phoenix, Arizona, USA          car  
 4 Little Rock, Arkansas, USA     car  
 5 Los Angeles, California, USA   car  
 6 San Francisco, California, USA car  
 7 Oakland, California, USA       car  
 8 Half Moon Bay, California, USA car  
 9 San Diego, California, USA     car  
10 Stamford, Connecticut, USA     car  
# ℹ 75 more rows
# ℹ Use `print(n = ...)` to see more rows

Plotting the us_travels dataset

For plotting the us_travels dataset use the following snippet can be used:

library(mapBliss)
data("us_travels")

locations <- us_travels$location
# The how argument needs to be one 
# less than the number of the locations
how <- us_travels$how[-1] 

plot_hybrid_route_flex(
    addresses = locations, 
    how = how,
    radius = 1,
    label_text = ""
  )

This yields the following leaflet map:

The map themes be visually altered by changing the mapBoxTemplate argument.

Conclusion

The us_travels dataset provides a nice example for how the {mapBliss} package can be used to make beautiful maps of your travels around the world from the R console with minimal effort.

For framing these maps, I presently recommend using a platform like Canva. I additionally recommend checking out the {starBliss} package for making custom star maps of the night sky.

Here is a visual I made for a print of a trip this past summer with {mapBliss}, {starBliss} and Canva:

There are some GitHub issues open related to making ready to use frames (which have had some success) but there’s still some considerable effort required to accomplish the desired goal. If you think you can do it, feel free to open a pull request!

Disclaimer: This was originally written on my Medium blog here, so the formatting is a little different from my usual style.

If you just got started or have been working a while in a data role, the jargon thrown around can sometimes get overwhelming with all the things to need to learn or be on top of. After spending some time thinking about it, I have managed to come up with 5 tips for landing a “data professional” role.

As martial arts skeptic Matt Thorton puts it when discussing the comparison of various martial arts:

“Just as there is no such thing as “Canadian geometry”, but there is geometry.”

While he wrote this in context of martial arts, I think the same is true with working with data. Business domains come and go, but many of the problems stay the same when working with data.

I am far from being a master in the fields of Data Science, Data Engineering and Machine Learning. However I have found the following useful in one setting or another in my experience.

With that in mind, my five tips are:

1. Learn to Pivot Data from Long to Wide Form (and vice versa)

This is something which comes up in analyzing or preparing data in one way or another. This is a skill which is language agnostic also. The more settings you know how to do it, the more versatile you will be in contributing to an analysis or data engineering project.

If you want to learn how to do this in R or Python, I have a blog on my personal site showing how to pivot data from long to wide from in R and Python and comparing the two.

Check out RvsPython #2: Pivoting Data From Long to Wide Form.

2. Learn to Manipulate Character Strings.

Being someone who specializes in data means being someone who is able to turn a mess into insights. Develop your regular expressions (regex) skills, learn to extract data from strings and blurbs of text.

Extracting data from strings is something that I have encountered regularly and I have done it in multiple roles. Its always gonna be needed. Get good at it.

There are a variety of ways which you can learn about regular expressions, but if you want something more structured I would recommend going through Hackerrank’s regex challenges.

Learning regular expressions might initially seem overwhelming with having to memorize everything, but I can tell you that you don’t need to memorize it all! Quite often I look at a string I want to extract data from and copy-paste it (or some of it) in regex101.com and figure out the regular expression I need through what I remember and some trial and error. The more you do it, the easier it gets!

3. Learn to Write API Requests and Handle the Response

This isn’t particularly hard, but oftentimes, data is not handed out on a silver platter for you to analyze or load. Mastering the art of data acquisition and transformation is required for all data professionals regardless of where they are in terms of workflow.

Get in the practice of reading API documentation. Create a request. See what the response is. 200- you’re in business. 502? 404? Not so much- figure out what the problem is and if there is a solution.

Its not incredibly difficult to learn how to write API requests and its not something which is language specific. Knowing how read API documentation and writing the appropriate requests is something I deal with regularly. I’m sure I’m not the only one out there. Learn how to do it.

For more info on understanding API requests DataCamp has a good article on the topic using the requests module in Python.

R user? Check out the CRAN documentation for getting started with thehttr package.

But what do you do if there is no accessible API available to get your data?

4. Learn to Scrape Data

Data rarely comes readily available. Often times you will need to figure out other ways to get it. For some reason web scraping is not spoken about enough. Its likely because there is moral hazard attached to it.

Provided that there is no TOS violations, web scraping enables collection of data to be optimized. Just be decent and not slam the servers with rapid speed requests. Put some system sleeps in between to have your code behave as a regular user as well.

Check out the selenium module. With this video it only takes around 10 minutes to learn and its a skill that will be worthwhile in a variety of settings.

With all of this in mind, lets get to the final point:

5. Remember: Machine Learning is useless without being able to control the data.

While having knowledge of machine learning is a powerful asset. If you are unable to manipulate and transform the data appropriately, the algorithms that you have are useless because they only accept data when shaped in a specific form.

If your company is looking to employ AI but does not have a data strategy, its likely that one needs to be set up before looking into using AI and machine learning.

Conclusion

The above are are some tips that are that I think are handy to be mindful of when being presented with recurring questions and problems relating to data as a data professional.

Be prepared!

I hope you found this blog useful. Be sure to share it around if you liked it and let me know if you didn’t!

I have been blogging for nearly two years and with it I got to learn a lot about writing good code, challenges that arise (particularly with data science) and how to navigate them.

What I also got to be exposed to the two truths below. Once I became aware of them, it helped a lot with filtering out the noise regarding exactly what I needed to do to be successful in this field. Its my hope that you will also.

The two truths are:

1. During any “gold rush” there will be people selling “picks and shovels”

During the California gold rush its said that those who made the most money were those who sold “picks and shovels”.

In today’s day and age the “gold rush” for many is a job in tech, and trying to get there comes with encountering many sellers of “picks and shovels”- be it in the form of a course, a software or a development framework.

This is not to say avoid paying for a course of a software that a individual or a company is selling. Do your own due diligence. There are many options for learning how to code that you can use for free that offer just a as good (if not better) learning experience than paid courses.

As far as learning is concerned, I am a huge fan of freeCodeCamp and the content it provides via full courses on Youtube.

If you ask me, in terms of educational content, freeCodeCamp is the one stop location for everything you need to know when started to learn a programming language or a specialization. I’ve even seen prep courses for networking and cloud computing certifications there!

If I could tell you one thing- before putting your money into sort of education or software, see if it does better than the free DIY option- you will be surprised how many aren’t!

2. Those who can’t do- Teach!

When most of us have heard this line it usually refers to insulting the profession of teaching. I’m using it for the exact opposite.

If you find yourself not competent in being able to solve a problem or write code without getting errors and spending hours of debugging- teach it!

Look into the problem and how to solve it. Share how you solved it. Share what you learned in the process in the form of a tweet, blog or video. Demonstrating what you can do speaks volumes more than a certificate of completion from many of the “courses” you’ll find for sale.

I am a big fan of having a blog which is not connected to a particular platform where you have complete control (like a WordPress site). That way you can migrate and customize your content at whenever you want.

Austin Kleon in his book “Show Your Work” writes about having a website which really reflected with me:

Don’t think of your website as a self-promotion machine, think of it as a self-invention machine.[…] Fill your website with your work and your ideas and the stuff you care about. Over the years, you will be tempted to abandon it for the newest, shiniest social network. Don’t give in. Don’t let it fall into neglect. Think about it in the long term. Stick with it, maintain it, and let it change with you over time.

If you aren’t experienced yet, having a blog sharing your progress will accelerate your ability to learn. Before you know it you will find people coming to you for help, and you will find yourself being able to answer their questions- and getting paid for it. Its happened for me and it can certainly happen for you.

Conclusion

These ideas came to me after spending more time and seeing people who were in my situation lost trying to figure out how to make their break in this field.

I hope you found this useful. Be sure to share this blog around!

Over the past while with my time on LinkedIn, I got to have exposure to many people from many different lines of work. I also managed to have carved a space for myself there where I can post about Data Science topics and share my blogs along the way. There have always been posts and polls comparing R and Python as well as the subsequent debates among users of the languages as far as which one is superior for doing Data Science. While these sort of arguments will never end and I am far from innocent of engaging in them, I chose to take to task understanding why Data Science practioners preferred one languge over another by “controlling” for exposure to the other language.

In this blog I am going to share my results from my LinkedIn polls comparing respondents preferences. The polls asked for respondents preferences for:

1. Using dplyr in R vs pandas in Python for data wrangling,

2. Using ggplot2 in R vs matplotlib and seaborn in Python for data visualization, and

3. Using Jupyter notebooks vs RMarkdown for writing reports.

Disclaimer

This is by no means a formal study, its more of just me sharing my findings in blog form. Social media platforms come and go, but having a blog where I can share my findings (albeit less popular) offers a place where I can post my curated content. Likely due to LinkedIn’s algorithms, my first and second questions got more traction with over 132,000 views combined and over 1600 and 1300 votes respectively, while my last question only got a little more than 4000 views and over 106 votes at the time of writing.

To quote a comment on one of my polls:

With this in mind, lets share the results of these polls.

(Visuals were made with ggplot2 and the ggtech package for the theme)

1. dplyr vs pandas

As expected, most users who were pro-pandas never used dplyr before. However, when controlling for prior experience, it was pretty much a 50-50 split among respondents between using pandas in Python and dplyr in R. There were some comments recommending that I check out the data.table and dtplyr packages in R; while I don’t have much exposure to using those packages presently, I hope to check them out in the future.

For my closest experience to dplyr in Python, check out my review on the siuba module.

2. ggplot2 vs matplotlib and seaborn

In the case of comparing ggplot2 to matplotlib and seaborn among users who had experience with both packages, ggplot2 is preferred by 56% of users. Most users of matplotlib and seaborn don’t have experience with ggplot2 and vice-versa.

I was told to check out the plotly library which is compatible in R and Python and it really looks like a great library to have for building interactive dashboards and applications. While I don’t have much experience with it now, I do hope to check it out when time allows for it.

3. Using Jupyter notebooks vs RMarkdown for writing reports.

The results from this poll are questionable as I only got 106 replies to this poll. With this in mind these are the results:

Of users with experience with using both RMarkdown and Jupyter notebooks for writing their reports, 63% of users prefer using RMarkdown over Jupyter notebooks, however there are more users who have experienced Jupyter notebooks than RMarkdown.

Conclusion

With all being said, using dplyr in R or pandas Python for doing data wrangling seems like a toss up among users with experience with both languages. For data visualization, ggplot2 seems to be preferred over matplotlib or seaborn and if you trust the sample size, RMarkdown is preferred over Jupyter notebook among users with experience with both.

In general, apparent that R is still the underdog in terms of it being a language used for Data Science and programming- but by no means does that make me intend on stopping from using it any time soon.

When I get the time, I look forward to giving data.table and plotly a spin!

Thank you for reading!

The following is the same content that I have posted in my other blog on this topic in R, but written in Python. While I actually first wrote the code for doing this in Python, I’ll be posting the similar verbiage from that blog here.

Introduction

Often when scraping data, websites will ask a user to enter a postal code to get the locations near it. If you are interested in collecting data on locations in Canada for an entire Province or the entire province from a site, it might be hard to find a list of all postal codes or FSAs in Canada or in a given province which is easy to use. Information on how FSAs work can be found here.

In this blog, I’m going to share a brief snippet of code that you can use to generate Canadian FSAs. While some FSAs generated may not actually exist, if we follow the rules about Canadian postal codes, it serves as a good substitute in lieu of an actual list.

The Code

is essentially 3 nested for-loops. While many programmers would not advise writing nested for-loops, I find that for this case it is easier to understand and write.

import string

fsa_list =[]

alphabet = list(string.ascii_uppercase)

province_alphabet = ["A","B","C","E","G","H","J","K","L","M","N","P","R","S","T","V","X","Y"]


nonLetters=["D", "F", "I", "O", "Q", "U"]

second_letter= list(set(string.ascii_uppercase)-set(nonLetters))

for letter1 in province_alphabet:
    for number in [0,1,2,3,4,5,6,7,8,9]:
        for letter2 in alphabet:
            fsa_list.append(letter1+str(number)+letter2)

We now have our simulated FSAs!

import random

random.sample(fsa_list,10)

['A0I', 'S7H', 'R5O', 'V4T', 'R9G', 'B6W', 'T6D', 'B9O', 'M0B', 'Y7H']

In this blog post, I am going to give the siuba module a spin by looking at the wine ratings dataset from the TidyTuesday project. While I hope to post an entire analysis in Python in the future, my focus is primarily on the features that the siuba has to offer. The questions, however trivial, are aimed at trying to highlight the features that are available.

The data

To stay in the spirit of doing things entirely in Python, I’m going to load the raw data from the github with pandas’ read.csv() function.

import pandas as pd


dt = pd.read_csv("https://raw.githubusercontent.com/rfordatascience/tidytuesday/master/data/2019/2019-05-28/winemag-data-130k-v2.csv")

Now that the data is loaded, lets begin with playing with the data! The following questions we’re going to ask are the following:

1. What are the 10 most frequent countries in the data set?
2. Who were the tasters who reviewed Cabernet Sauvignon wines?
3. What was is the 5 number summary for cost of Merlot wines?
4. What were the top wine varieties with the best scores?

While these questions are basic, it will give us the opportunity to show what the siuba module can do. I’m sure there are more sophisticated questions that we can ask, but this should be good for now to show off what can be done with siuba.

So without further ado, lets get into it!

1. The 10 most frequent countries in the data set.

This seems like something pretty standard that can be done in R with dplyr. The following is a comparison between dplyr in R and siuba in Python:

With dplyr (in R)

dt2 %>% 
  select(country) %>% 
  group_by(country) %>%
  summarize(n=n()) %>% 
  arrange(desc(n)) %>% 
  top_n(10)

## Error: `n()` must only be used inside dplyr verbs.

With siuba (in Python)

from siuba import *
from siuba.dply.vector import *

(
  dt >> 
  select(_.country)>>
  group_by(_.country)>>
  summarize(n= n( _.country)) >>

  # using desc
  arrange(desc(_.n)) >>
  top_n(10)
)

     country      n
40         US  54504
37      Spain   6645
31   Portugal   5691
22      Italy  19540
17    Germany   2165
15     France  22093
8       Chile   4472
3     Austria   3345
2   Australia   2329
0   Argentina   3800

There are a few differences that are immediately noticeable. With siuba we have to import all the the functions available and also import the siuba.dply.vector to get the relevant dplyr functions. Additionally, the desc(_.n) argument in the arrange() doesn’t give the top 10 most frequented countries in descending order.

After a bit of trial and error and some Googling, to remedy the siuba code, instead of writing arrange(desc(_.n)), we can write arrange(-_.n).

(
dt >> 
  select(_.country)>>
  group_by(_.country)>>
  summarize(n= n( _.country)) >>

  # desc doesn't work
  arrange(-_.n) >>
  top_n(10)

)

      country      n
40         US  54504
15     France  22093
22      Italy  19540
37      Spain   6645
31   Portugal   5691
8       Chile   4472
0   Argentina   3800
3     Austria   3345
2   Australia   2329
17    Germany   2165

Much better!

2.The tasters who reviewed Cabernet Sauvignon wines

The difference between siuba and dplyr here seems to be more rooted in the difference between R and Python (but I could be wrong). R accounts for NA values while Python does not.

With dplyr (in R)

dt2 %>% 
  filter(variety =="Cabernet Sauvignon") %>% 
  select(taster_name) %>% 
  group_by(taster_name) %>% 
  summarize(n=n()) %>% 
  arrange(desc(n))

## Error: `n()` must only be used inside dplyr verbs.

With siuba (in Python)

(
  dt >>
  filter(_.variety=="Cabernet Sauvignon") >>
  select(_.taster_name)>>
  group_by(_.taster_name)>>
  summarize(n= n( _.taster_name)) >>

  # using desc
  arrange(-_.n)
)

           taster_name     n
17      Virginie Boone  1829
11   Michael Schachner  1358
15    Sean P. Sullivan   785
13        Paul Gregutt   677
10       Matt Kettmann   454
6           Jim Gordon   399
7       Joe Czerwinski   263
9        Lauren Buzzeo   140
1   Anna Lee C. Iijima    94
16     Susan Kostrzewa    67
12       Mike DeSimone    66
14          Roger Voss    49
5         Jeff Jenssen    35
8        Kerin O'Keefe    32
0   Alexander Peartree     8
3         Carrie Dykes     7
4          Fiona Adams     6
2     Anne Krebiehl MW     2

What I like about R is that it accounts for missing data which I can see and account for, but with the way I wrote the Python code I am unable to see this.

3. 5 number summary for price of Merlot wines

With dplyr (in R)

dt2 %>% 
  filter(variety =="Merlot") %>% 
  select(price) %>% 
  summary()

##      price       
##  Min.   :  4.00  
##  1st Qu.: 15.00  
##  Median : 24.00  
##  Mean   : 29.54  
##  3rd Qu.: 35.00  
##  Max.   :625.00  
##  NA's   :40

With siuba (in Python)

(
dt >> 
  filter(_.variety =="Merlot")>>
  select(_.price)
).describe()

             price
count  3062.000000
mean     29.543436
std      33.340882
min       4.000000
25%      15.000000
50%      24.000000
75%      35.000000
max     625.000000

With dplyr I am able to pipe my filtered data directly into summary() to get a five number summary. With siuba, since the object output is still a pandas object (for good reason), I am only able to get the 5 number summary by using the .describe() method.

4. Top wine varieties with the best scores

With dplyr (in R)

dt2 %>% 
  distinct(variety,points) %>% 
  filter(points==max(points))

## # A tibble: 13 x 2
##    points variety                   
##     <dbl> <chr>                     
##  1    100 Muscat                    
##  2    100 Prugnolo Gentile          
##  3    100 Champagne Blend           
##  4    100 Merlot                    
##  5    100 Portuguese Red            
##  6    100 Sangiovese                
##  7    100 Cabernet Sauvignon        
##  8    100 Bordeaux-style Red Blend  
##  9    100 Chardonnay                
## 10    100 Sangiovese Grosso         
## 11    100 Syrah                     
## 12    100 Port                      
## 13    100 Bordeaux-style White Blend

With siuba (in Python)

import numpy as np
(
dt >> 
  distinct(_.variety,_.points) >>
  filter(_.points==np.max(_.points))
)

                         variety  points
186                       Muscat     100
1188            Prugnolo Gentile     100
2477             Champagne Blend     100
2539                      Merlot     100
2629              Portuguese Red     100
2744                  Sangiovese     100
2747          Cabernet Sauvignon     100
3047    Bordeaux-style Red Blend     100
3567                  Chardonnay     100
3840           Sangiovese Grosso     100
3875                       Syrah     100
3887                        Port     100
3977  Bordeaux-style White Blend     100

Aside for numpy being imported to get the max() function the code with siuba and dplyr is the same.

Conclusion

With everything I checked out. Aside for a hiccup with the desc() function, using the siuba module is great for making use of dplyr verbs with dataframes in Python. If properly maintained and stable, I can see tons of applications of siuba being used beyond jupyter notebooks in back end programs which deal with data.

This is only my small review, so I recommend everyone else give it a try themselves by installing it (pip install siuba) and seeing what you think! You won’t be disappointed!

Because I thought it was really cool to see live prices of stocks/crypto in real time on Discord, I wanted to make one myself. Lets run through the steps for making a ticker bot like I did today!

The script requires that you have the discord and the yahoo_fin modules installed. So before you start be sure to have them installed in your environment.

pip install discord
pip install yahoo_fin

Step 1: Go to the Discord Developers Portal and create an application with a bot

First things first, you’re going to want to log into the Discord developers portal navigate to “Applications” and click on “New Application” after which you will be given a prompt to put in the name of the application.

After the application is created, click on your application and create a bot by going to the “Bot” section of the settings and create a bot.

After that, go to the OAuth2 section of the bot and under “scopes” click “bot” and invite the bot with the link generated.

Once you invited you can proceed with developing your own bot.

Step 2: Get your token

It is very important to have your unique token for your bot. This can be found in your application settings under the “bot” tab.

Once the token is copied, put it in a .txt file that you can access.

Step 3: The Script

The script follows the following structure:

from discord.ext import commands
from yahoo_fin import stock_info as si


def read_token():
    with open('botTokenPath.txt', 'r') as f:
        lines = f.readlines()
        return lines[0].strip()

token = read_token()

bot = commands.Bot(command_prefix="!!")

@bot.command()

async def run(ctx):
    while True:
        change = si.get_quote_data("StockSymbol")

        name="S&P500: \n$"+ str(round(si.get_live_price("StockSymbol"),2))+ "("+str(round(change["regularMarketChangePercent"],2))+")"+"%"
        await ctx.guild.me.edit(nick=name)

Run the script, and there you have it! Once you have your token and Stock Symbols picked, type !run in one of your channel servers and you should have your stock ticker up and running!

As for hosting, there are lots of great resources out there for figuring how to do that. If you’re curious, I recommend either the solution offered by FreeCodeCamp (see the end of the article) or by Tech With Tim.

Safety Tips

1. Never share your bot token with anyone– this gives the opportunity for hackers to enter into your server and expose vulnerabilities and ruin your overall discord experience
2. Be care with the permissions you give your bot– the permissions you give your bot will enable your bot to function. But if you give any bot full reign to do whatever it wants-regardless of the abilities you programmed it to have- you are exposing your server to being hacked or nuked . For these bots, the only real permission they should have is the ability to change their nicknames. Everything else is not needed.

Thank you for reading! If you have any questions be sure to reach out!

For the technical background I took the introductory material from the paper “Convolutional neural networks: an overview and application in radiology” by Yamashita, Nishio, Do and Togashi from the Journal “Insights into Imaging” (Published June 27, 2018). For the applications I produced some original analysis using some of the toy data sets provided by keras (MNIST numbers and the CFAIR-10 datasets) and compared how convolutional neural networks preformed relative to simple feed forward networks.

Concluding remarks were based on the publicly available jupyter notebook written by well known python-programming Youtuber Tech with Tim here (if you are reading this – thank you for the resource Tim!)

As far as the actual Jupyter Notebook with the code I used, when I get some time I hope to reformat it and have it available to see!

Below are the slides I used for this talk. Let me know what you think!

In this blog I’m going to do the following:

1. Perform Monte Carlo Using R and Python with for loops.
2. Use some of Python’s best practices and see how it compares to R’s lapply.

(Note: There is already a popular article article comparing for loops in R and Python against R’s lapply here)

Disclaimer

This is a follow up blog to my last write up comparing R and Python with Monte Carlo simulations. For context, check out the blog here.

Using for loops in R

In the last post, we didn’t have our generated data timed. I thought it would be a good idea to include it in the total processing time. I will also be changing the method I generated the data by using R’s for loop and will start timing the code from there.

I did my best to make my R code as similar to the Python code in the last blog- if you see an issue, please comment!

#' Define Number of points we want to estimate
n<-c(10,100,1000,10000,100000,1000000)

#' Our Transformation function

y<- function(u) {
  4*sqrt(1-u^2)
}

#' Start the timer

startTime<-Sys.time()

#' Generate our random uniform variables
x<-list()

for(i in 1:length(n)){
  x[[i]]<-runif(n[i])
}

#' Transform our uniform variables.

yvals<-list()
for (i in 1:length(x)){
  #' Need to define this so that the list element will be populated
  #' See: https://stackoverflow.com/questions/14333525/error-in-tmpk-subscript-out-of-bounds-in-r
  yvals[[i]]<-1
  for(j in 1:length(x[[i]])){
    yvals[[i]][j]<-y(x[[i]][j])
  }
}


#' Calculate our approximations of pi

avgs<- c()

for(i in 1:length(yvals)){
  avgs[i]<-mean(yvals[[i]])
}

endTime<-Sys.time()-startTime

endTime

## Time difference of 1.009413958 secs

data.frame(n, "MC Estimate"=unlist(avgs), "Difference from True Pi"= abs(unlist(avgs)-pi))

##         n MC.Estimate Difference.from.True.Pi
## 1      10 3.281637132         0.1400444782036
## 2     100 3.391190973         0.2495983193740
## 3    1000 3.090265904         0.0513267494211
## 4   10000 3.143465663         0.0018730098616
## 5  100000 3.141027069         0.0005655842822
## 6 1000000 3.141768899         0.0001762457079

Using for loops in Python (From previous blog)

As I did in the previous blog, here is the code I used to run the Monte Carlo algorithm with for loops. I heard there are more accurate ways to time this code, but since I want it to be similar to my R code- I am doing it this way.

import numpy as np
import pandas as pd
import time

# Define Number of points we want to estimate

n = [10, 100, 1000, 10000, 100000, 1000000]

# Our Transformation function

def y(x):
    return 4 * np.sqrt(1 - x ** 2)


#Start the timer

startTime= time.time()

# Generate our random uniform variables

x = [np.random.uniform(size=n) for n in n]



startTime= time.time()
yvals = []
for array in x:
    yval=[]
    for i in array:
        yval.append(y(i))
    yvals.append(yval)

avgs=[]

for array in yvals:
  avgs.append(np.mean(array))

endTime= time.time()-startTime

# How long it took to run our code
print("Time difference of "+ str(endTime) + " secs\n")


# Output

## Time difference of 2.790393352508545 secs

print("Estimated Values of Pi\n")

## Estimated Values of Pi

pd.DataFrame({"n":n,
              "MC Estimate":avgs,
              "Difference from True Pi": [np.abs(avg-np.pi) for avg in avgs]})

##          n  MC Estimate  Difference from True Pi
## 0       10     3.259405                 0.117812
## 1      100     3.351556                 0.209963
## 2     1000     3.130583                 0.011009
## 3    10000     3.126542                 0.015050
## 4   100000     3.144484                 0.002891
## 5  1000000     3.140740                 0.000853

library(reticulate)

py$endTime/as.numeric(endTime)

## [1] 2.764369693

Ok- so using for loops R isn’t as fast as I initally stated. However, based on my machine R is still over twice as fast as Python with for loops.

Hey, it ain’t 220 but its something

Using R’s “Best Practices” (Using the apply family)

Instead of using for loops, a faster alternative is to use the apply family of functions, namely sapply and lapply.

#' Start the timer
startTime<-Sys.time()

#' Generate our random uniform variables
x<-sapply(n,runif)

yvals<-lapply(x,y)
avgs<-lapply(yvals,mean)

newendTime<-Sys.time()-startTime
newendTime

## Time difference of 0.1879060268 secs

#' Speed - for loop vs apply
as.numeric(endTime)/as.numeric(newendTime)

## [1] 5.371908366

Using Python’s best practices

After getting several comments of (constructive) criticism about how the comparison was not fair here’s some new code implementing some of the best practices in writing faster code.

This is some code that I saw posted in a comment on my LinkedIn post (thank you Thomas Halvorson), which is pretty similar in structure to the R code I have listed above.

I’m sure there are better ways out there (I have seen in the comments for the last blog a lot of very good solutions), but I found this to be the most readable and follows a structure similar to R’s.

(Let me know if you have something better!)

startTime= time.time()

x = [np.random.uniform(size=n) for n in n]
yvals = list(map(y, x))
avgs = list(map(np.mean, yvals))

endTime= time.time()-startTime

# How long it took to run our code
print("Time difference of "+ str(endTime) + " secs\n")

## Time difference of 0.0629582405090332 secs

Comparing R with Python now we have:

as.numeric(newendTime)/py$endTime

## [1] 2.984613695

Python is nearly 3 times faster on my machine using the updated code.

Conclusion

Well, you live and learn. Best practices can make it or break it for your code and this updated analysis can help give you a better idea.

Thank you everyone for reading my last blog post and pointing out some obvious issues that I didn’t notice! I definitely will be using the map() function more often in my Python code!

Did you like this content? Be sure to never miss an update and Subscribe!

This blog has got alot of criticism as to the lack of best practices being used in the Python code. I have written a follow up to this issue here. Let me know what you think!

Introduction

The number , while being an irrational and transcendental number is a central number to which much of mathematics and science at large relies on for many calculations.

A whole book can be written on this matter alone, but today we are going to focus on approximating the true value of using Monte Carlo simulations in R and Python!

Disclaimer

This problem has been covered extensively across the internet and serves as a benchmark example of what Monte Carlo can do. What we are going to do is highlight how this method preforms in both R and Python.

The General Algorithm

The formula for the unit circle is:

For , the length of an eighth of the unit circle is . Thus, to approximate , the function that we will be using is:

To approximate , we use the following Algorithm:

1. Generate
2. Calculate
3. Take the mean of , : – this is our approximation

The code for this is relatively straight forward. But the question is, which code will run faster?

Let’s go!

The Test

For our challenge we are going to be writing code which is as intuitive as possible in each language. We are going to seek to approximate the value of using the above algorithm. For R we will be using lapply to implement the Monte Carlo algorithm and for Python we will using for loops.

The Solution with R

#' Define Number of points we want to estimate
n<-c(10,100,1000,10000,100000,1000000)

#' Generate our random uniform variables
x<-sapply(n,runif)

#' Our Transformation function

y<- function(u) {
  4*sqrt(1-u^2)
}


startTime<-Sys.time()
yvals<-lapply(x,y)
endTime<-Sys.time()-startTime
avgs<-lapply(yvals,mean)

endTime

## Time difference of 0.01399588585 secs

data.frame(n, "MC Estimate"=unlist(avgs), "Difference from True Pi"= abs(unlist(avgs)-pi))

##         n MC.Estimate Difference.from.True.Pi
## 1      10 3.281637132         0.1400444782036
## 2     100 3.391190973         0.2495983193740
## 3    1000 3.090265904         0.0513267494211
## 4   10000 3.143465663         0.0018730098616
## 5  100000 3.141027069         0.0005655842822
## 6 1000000 3.141768899         0.0001762457079

The Solution with Python

import numpy as np
import pandas as pd
import time

# Define Number of points we want to estimate

n = [10, 100, 1000, 10000, 100000, 1000000]

# Generate our random uniform variables

x = [np.random.uniform(size=n) for n in n]


# Our Transformation function

def y(x):
    return 4 * np.sqrt(1 - x ** 2)


startTime= time.time()
yvals = []
for array in x:
    yval=[]
    for i in array:
        yval.append(y(i))
    yvals.append(yval)

avgs=[]

for array in yvals:
  avgs.append(np.mean(array))

endTime= time.time()-startTime

# How long it took to run our code
print("Time difference of "+ str(endTime) + " secs\n")


# Output

## Time difference of 3.146182060241699 secs

## Estimated Values of Pi

pd.DataFrame({"n":n,
              "MC Estimate":avgs,
              "Difference from True Pi": [np.abs(avg-np.pi) for avg in avgs]})

##          n  MC Estimate  Difference from True Pi
## 0       10     3.320525                 0.178933
## 1      100     3.172290                 0.030698
## 2     1000     3.156044                 0.014451
## 3    10000     3.141675                 0.000083
## 4   100000     3.147255                 0.005662
## 5  1000000     3.141400                 0.000193

Comparing R with Python

From the following ratio we can see how much faster R is than Python:

library(reticulate)

reticulate::py$endTime/as.numeric(endTime)

## [1] 224.7933496

Woah! Using my machine- R is over 220 times faster than Python!

I think it’s pretty clear to see who the winner is as far as speed is concerned.

Concluding Remarks

While R most of the time sits on the sidelines in the Python-dominant world of Data Science- we need to keep in mind where Python’s weaknesses lie and when to pivot from and use R.

Doing simulation? R, please.

Did you like this content? Be sure to never miss an update and Subscribe!