(c) 2019 Justin Bois. With the exception of pasted graphics, where the source is noted, this work is licensed under a Creative Commons Attribution License CC-BY 4.0. All code contained herein is licensed under an MIT license.
This document was prepared at Caltech with financial support from the Donna and Benjamin M. Rosen Bioengineering Center.
This lesson was generated from a Jupyter notebook. You can download the notebook here.
import os
import glob
Reading data in from files and then writing your results out again is one of the most common practices in scientific computing. In this tutorial, we will learn about some of Python's File I/O capabilities. We will use a PDB file as an example. The PDB file contains the crystal structure for the tetramerization domain of p53.It is stored in the file ~/git/bootcamp/data/1OLG.pdb
. (Make sure you launch your notebook from the ~/git/bootcamp/
directory.) Note that 1OLG
is its unique Protein Databank identifier.
To open a file, we use the built-in open()
function. When opening files, we should do this using context management. I will demonstrate how to open a file and then describe the syntax.
with open('data/1OLG.pdb', 'r') as f:
print(type(f))
Python has a wonderful keyword, with
. This keyword enables context management. Upon entry into a with
block, variables have certain meaning. In this case, the variable f
has the meaning of an open file, an instance of the _io.TextIOWrapper
class. Upon exit, certain operations take place. For file objects created by opening them, the file is automatically closed upon exit, even if there is an error. This is important. If your program raises an exception before you have a chance to close the file, it won't get closed and you could be in trouble. If you use context management, the file will still get closed. So here is an important tip:
Let's focus for a moment on the variable f
in the above code cell. It is a Python file
object, which has methods and attributes, just like any other object. We'll explore those in a moment, but first, let's look at how we opened the file. The first argument to open()
is a string that has the name of the file, with the full path if necessary. The second argument is a string that says what we will be doing with the file. I.e., are we reading or writing to the file? The possible strings for this second argument are
string | meaning |
---|---|
'r' |
open a text file for reading |
'w' |
create and open a text file for writing |
'a' |
append an existing text file |
'r+' |
open a text file for reading and writing |
append 'b' to any of the above |
same as above, except for binary files |
We will mostly be working with text files in the bootcamp, so the first three are the most useful. A big warning, though....
We will focus on the methods f.read()
and f.readlines()
. What do they do?
method | task |
---|---|
f.read() |
Read the entire contents of the file into a string |
f.readlines() |
Read the entire file into a list with each item being a string representing a line |
First, we'll try using the first method to get a single string with the entire contents of the file.
# Read file into string
# Read contents of the file in as a list
with open('data/1OLG.pdb', 'r') as f:
f_str = f.read()
# Let's look at the first 1000 characters
f_str[:1000]
We see lots of \n
, which signifies a new line. The backslash is known as an escape character, meaning that the n
after it does not signify the letter n, but that \n
together means a new line.
Now, let's try reading it in as a list.
# Read contents of the file in as a list
with open('data/1OLG.pdb', 'r') as f:
f_list = f.readlines()
# Look at the list (first ten entries)
f_list[:10]
We see that each entry is a line, including the newline character. To look at lines in files, the rstrip()
method for strings can come it handy. It strips all whitespace, including newlines, from the end of a string.
f_list[0].rstrip()
What if we do not want to read the entire file into a list? For example, if a file is several gigabytes, we do not want to spend all of our RAM storing a list. Instead, we can read it line-by-line. Conveniently, the file object can be used as an iterator.
# Print the first ten lines of the file
with open('data/1OLG.pdb', 'r') as f:
for i, line in enumerate(f):
print(line.rstrip())
if i >= 10:
break
Alternatively, we can use the method f.readline()
to read a single line in the file and return it as a string.
# Print the first ten lines of the file
with open('data/1OLG.pdb', 'r') as f:
i = 0
while i < 10:
print(f.readline().rstrip())
i += 1
Each subsequent call to f.readline()
reads in the next line of the file. (Remember, as we read through a file, we keep moving forward in the bytes of the file and we have to use f.seek()
to rewind.)
Writing to a file has similar syntax. We already saw how to open a file for writing. Again, context management is useful. However, before trying to open a file, we should check to make sure a file of the same name does not exist before opening it. The os.path
module is useful. The function os.path.isfile()
function checks to see if a file exists.
os.path.isfile('data/1OLG.pdb')
Now that we know how to check existence of a file so we do not overwrite it, we can open and write a file.
if os.path.isfile('mastery.txt'):
raise RuntimeError('File mastery.txt already exists.')
with open('mastery.txt', 'w') as f:
f.write('This is my file.')
f.write('There are many like it, but this one is mine.')
f.write('I must master my file like I must master my life.')
Note that we can use the f.write()
method to write strings to a file. Let's look at the file contents.
!cat mastery.txt
Ah! There are no newlines! When writing to a file, unlike when you use the print()
function, you must include the newline characters. Let's try again, intentionally obliterating our first attempt.
with open('mastery.txt', 'w') as f:
f.write('This is my file.\n')
f.write('There are many like it, but this one is mine.\n')
f.write('I must master my file like I must master my life.\n')
!cat mastery.txt
That's better. Note also that f.write()
only takes strings as arguments. You cannot pass numbers. They must be converted to strings first.
# This will result in an exception
with open('gimme_phi.txt', 'w') as f:
f.write('The golden ratio is φ = ')
f.write(1.61803398875)
Yup. It must be a string. Let's try again.
with open('gimme_phi.txt', 'w') as f:
f.write('The golden ratio is φ = ')
f.write('{phi:.8f}'.format(phi=1.61803398875))
!cat gimme_phi.txt
That works!
As an example on how to do file I/O, we will take the PDB file and extract only the ATOM
records for the first chain of the tetramer and write only those entries to a new file.
It is useful to know that according to the PDB format specification, column 21 in the ATOM
entry gives the ID of the chain.
We also conveniently use the fact that we can have multiple files open in our with
block, separating them with commas.
with open('data/1OLG.pdb', 'r') as f, open('atoms_chain_A.txt', 'w') as f_out:
# Put the ATOM lines from chain A in new file
for line in f:
if len(line) > 21 and line[:4] == 'ATOM' and line[21] == 'A':
f_out.write(line)
Let's see how we did!
!head -10 atoms_chain_A.txt
!tail -10 atoms_chain_A.txt
Nice!
In the above snippet of code, we extracted all atom records from a PDB file. We might want to do this (or some other operation) for many files. For example, the directory ~/git/data/
has four PDB files in it. For the present discussion, let's say we want to pull the sequence of chain A out of each PDB file.
The glob
module from the standard library enables us to get a list of all files that match a pattern. In our case, we want all files matching data/*.pdb
, where *
is a wild card character, meaning that any matches of characters where *
appears are allowed. Let's see what glob.glob()
gives us.
file_list = glob.glob('data/*.pdb')
file_list
We have the four PDB files. We can now loop over them and pull out the sequences.
# Dictionary to hold sequences
seqs = {}
# Loop through all matching files
for file_name in file_list:
# Extract PDB ID
pdb_id = file_name[file_name.find('/')+1:file_name.rfind('.')]
# Initialize sequence string, which we build as we go along
seq = ''
with open(file_name, 'r') as f:
for line in f:
if len(line) > 11 and line[:6] == 'SEQRES' and line[11] == 'A':
seq += line[19:].rstrip() + ' '
# Build sequence with dash-joined three letter codes
seq = '-'.join(seq.split())
# Store in the dictionary
seqs[pdb_id] = seq
Let's take a look at what we got. We'll look at actin.
seqs['1J6Z']
%load_ext watermark
%watermark -v -p jupyterlab