This page illustrates how to use Python to perform a simple but complete analysis: retrieve data, do some computations based on it, and visualise the results.
Don’t worry if you don’t understand everything on this page! Its purpose is to give you an example of things you can do and how to go about doing them - you are not expected to be able to reproduce an analysis like this in Python at this stage! We will be looking at the concepts and practices introduced on this page as we go along the course.
As we show the code for different parts of the work, we will be touching on various aspects you may want to keep in mind, either related to Python specifically, or to research programming more generally.
We can use programs for our entire research pipeline. Not just big scientific simulation codes, but also the small scripts which we use to tidy up data and produce plots. This should be code, so that the whole research pipeline is recorded for reproducibility. Data manipulation in spreadsheets is much harder to share or check.
You can see another similar demonstration on the software carpentry site. We’ll try to give links to other sources of Python training along the way. Part of our approach is that we assume you know how to use the internet! If you find something confusing out there, please bring it along to a RITS drop-in session(UCL colleagues only).
In this course, we’ll always try to draw your attention to other sources of information about what we’re learning. Paying attention to as many of these as you need to, is just as important as these core notes.
Research programming is all about using libraries: tools other people have provided that do many cool things. By combining them we can feel really powerful but doing minimum work ourselves. The python syntax to import someone else’s library is “import”.
import geopy # A python library for investigating geographic information. # https://pypi.python.org/pypi/geopy~
Now, if you try to follow along on this example in an IPython notebook, you’ll probably find that you just got an error message.
You’ll need to wait until we’ve covered installation of additional python libraries later in the course, then come back to this and try again. For now, just follow along and try get the feel for how programming for data-focused research works.
geocoder = geopy.geocoders.Yandex(lang="en_US") geocoder.geocode('Cambridge', exactly_one=False)
[Location(Cambridge, Cambridgeshire County, United Kingdom, (52.208145, 0.133023, 0.0)),
Location(Cambridge, Ontario, Canada, (43.370599, -80.318989, 0.0)),
Location(Cambridge, Middlesex County, Massachusetts, United States of America, (42.385899, -71.141684, 0.0)),
Location(Cambridge, Washington County, State of Idaho, United States of America, (44.572039, -116.678507, 0.0)),
Location(Cambridge, Isanti County, Minnesota, United States of America, (45.572595, -93.223783, 0.0)),
Location(Cambridge, East London, Republic of South Africa, (-32.978267, 27.884695, 0.0)),
Location(Cambridge, Waikato, New Zealand, (-37.889307, 175.465014, 0.0)),
Location(Cambridge, Cambridgeshire County, United Kingdom, (52.210303, 0.176447, 0.0)),
Location(Cambridge Bay, Nunavut, Canada, (69.11673, -105.067782, 0.0)),
Location(Cambridge, Saint James, Jamaica, (18.291699, -77.895883, 0.0))]
The results come out as a list inside a list: [Name, [Latitude, Longitude]]. Programs represent data in a variety of different containers like this
Code after a # symbol doesn’t get run.
print("This runs") # print ("This doesn't") # print ("This doesn't either")
We can wrap code up in a function, so that we can repeatedly get just the information we want.
def geolocate(place): return geocoder.geocode(place, exactly_one = False)
Defining functions which put together code to make a more complex task seem simple from the outside is the most important thing in programming. The output of the function is stated by “return”; the input comes in brackets after the function name:
We can store a result in a variable:
london_location = geolocate("London") print(london_location)
The Yandex API allows us to fetch a map of a place, given a longitude and latitude. The URLs look like: https://static-maps.yandex.ru/1.x/?size=400,400&ll=-0.1275,51.51&z=10&l=sat&lang=en_US We’ll probably end up working out these URLs quite a bit. So we’ll make ourselves another function to build up a URL given our parameters.
import requests def request_map_at(lat, long, satellite=True, zoom=10, size=(400, 400)): base = "https://static-maps.yandex.ru/1.x/?" params = dict( z = zoom, size = str(size) + "," + str(size), ll = str(long) + "," + str(lat), l = "sat" if satellite else "map", lang = "en_US" ) return requests.get(base,params=params)
map_response = request_map_at(51.5072, -0.1275)
Let’s see what URL we ended up with:
url = map_response.url print(url[0:50]) print(url[50:100]) print(url[100:])
We can write automated tests so that if we change our code later, we can check the results are still valid.
from nose.tools import assert_in assert_in("https://static-maps.yandex.ru/1.x/?", url) assert_in("ll=-0.1275%2C51.5072", url) assert_in("z=10", url) assert_in("size=400%2C400", url)
Our previous function comes back with an Object representing the web request. In object oriented programming, we use the .
operator to get access to a particular property of the object, in this case, the actual image at that URL is in the
content property. It’s a big file, so I’ll just get the first few chars:
I’ll need to do this a lot, so I’ll wrap up our previous function in another function, to save on typing.
def map_at(*args, **kwargs): return request_map_at(*args, **kwargs).content
I can use a library that comes with Jupyter notebook to display the image. Being able to work with variables which contain images, or documents, or any other weird kind of data, just as easily as we can with numbers or letters, is one of the really powerful things about modern programming languages like Python.
import IPython map_png = map_at(*london_location)
print("The type of our map result is actually a: ", type(map_png))
The type of our map result is actually a: <class ‘bytes’>
Now we get to our research project: we want to find out how urbanised the world is, based on satellite imagery, along a line between two cites. We expect the satellite image to be greener in the countryside.
We’ll use lots more libraries to count how much green there is in an image.
from io import BytesIO # An object to convert between files and unicode strings import numpy as np # A library to deal with matrices import imageio # A library to deal with images
Let’s define what we count as green:
def is_green(pixels): threshold=1.1 greener_than_red = pixels[:,:,1] > threshold* pixels[:,:,0] greener_than_blue = pixels[:,:,1] > threshold*pixels[:,:,2] green = np.logical_and(greener_than_red, greener_than_blue) return green
This code has assumed we have our pixel data for the image as a 400 x 400 x 3 3-d matrix, with each of the three layers being red, green, and blue pixels.
We find out which pixels are green by comparing, element-by-element, the middle (green, number 1) layer to the top (red, zero) and bottom (blue, 2)
Now we just need to parse in our data, which is a PNG image, and turn it into our matrix format:
def count_green_in_png(data): f = BytesIO(data) pixels = imageio.imread(f) # Get our PNG image as a numpy array return np.sum(is_green(pixels))
We’ll also need a function to get an evenly spaced set of places between two endpoints:
def location_sequence(start, end, steps): lats = np.linspace(start, end, steps) # "Linearly spaced" data longs = np.linspace(start,end, steps) return np.vstack([lats, longs]).transpose()
location_sequence(geolocate("London"), geolocate("Cambridge"), 5)
array([[ 5.15073510e+01, -1.27660000e-01],
[ 5.16825495e+01, -6.24892500e-02],
[ 5.18577480e+01, 2.68150000e-03],
[ 5.20329465e+01, 6.78522500e-02],
[ 5.22081450e+01, 1.33023000e-01]])
We should display the green content to check our work:
def show_green_in_png(data): pixels = imageio.imread(BytesIO(data)) # Get our PNG image as rows of pixels green = is_green(pixels) out = green[:, :, np.newaxis] * np.array([0, 1, 0])[np.newaxis, np.newaxis, :] buffer = BytesIO() result = imageio.imwrite(buffer, out, format='png') return buffer.getvalue()
IPython.core.display.Image( map_at(*london_location, satellite=True) )
IPython.core.display.Image( show_green_in_png( map_at( *london_location, satellite=True)))
We can loop over each element in out list of coordinates, and get a map for that place:
for location in location_sequence(geolocate("London"), geolocate("Birmingham"), 4): IPython.core.display.display( IPython.core.display.Image(map_at(*location)))
So now we can count the green from London to Birmingham!
[count_green_in_png(map_at(*location)) for location in location_sequence(geolocate("London"), geolocate("Birmingham"), 10)]
[81911, 120660, 129892, 116438, 72826, 52788, 56890, 56208, 73858, 99759]
Let’s plot a graph.
import matplotlib.pyplot as plt %matplotlib inline
plt.plot([count_green_in_png(map_at(*location)) for location in location_sequence(geolocate("London"), geolocate("Birmingham"), 10)])
[<matplotlib.lines.Line2D at 0x2b7862b3f0f0>]
From a research perspective, of course, this code needs a lot of work. But I hope the power of using programming is clear.
We built little pieces of useful code, to:
By putting these together, we can make a function which can plot this graph automatically for any two places:
def green_between(start, end,steps): return [count_green_in_png(map_at(*location)) for location in location_sequence( geolocate(start), geolocate(end), steps)]
plt.plot(green_between('New York', 'Chicago', 20))
« [<matplotlib.lines.Line2D at 0x2b7862bde278>]
And that’s it! We’ve covered, very very quickly, the majority of the python language, and much of the theory of software engineering.
Now we’ll go back, carefully, through all the concepts we touched on, and learn how to use them properly ourselves.