Exercise¶
assignment 1: similarity of sequences
In this exercise you will write a program that calculates the distance between two sequences, e.g:
A simple program (without fonctions or modules) is sufficient.Solution
# define sequences
seq1 = "ACGT"
seq2 = "AGGT"
# initiate counter
distance_score = 0
# for each letter in the sequences
for a,b in zip(seq1, seq2)
# if they don't match, add a distance point
if a != b:
distance_score += 1
# print result to the terminal
print("Distance between A and B: ", distance_score)
-
Calculate the distance between the following sequences and print out the result. Since the following sequences are already aligned, we can calculate the distance between them. Change your program so that it can read two aligned sequences from the command line. Test your program with the following sequences.
Solution
import sys
# read sequences from command line arguments
seq1 = sys.argv[1]
seq2 = sys.argv[2]
# initiate counter
distance_score = 0
# for each letter in the sequences
for a,b in zip(seq1, seq2)
# if they don't match, add a distance point
if a != b:
distance_score += 1
# print result to the terminal
print(f"Distance between A and B: {distance_score}")
- Extend the program that the aligned sequences are printed out additionally to their distance.
Solution
-
Extend the program that the distance between two sequences is only calculated when both sequences have the same length. Test your program with the input sequences:
Note: Here you can use either the
assert
keyword or anif
clause. You could also raise a custom exception usingtry-except
.
Solution
import sys
# read sequences from command line arguments
seq1 = sys.argv[1]
seq2 = sys.argv[2]
assert len(seq1) == len(seq2), "Sequences must be of same length"
# initiate counter
distance_score = 0
# for each letter in the sequences
for a,b in zip(seq1, seq2)
# if they don't match, add a distance point
if a != b:
distance_score += 1
# print result to the terminal
print(f"Distance between A and B: {distance_score}")
-
Extend the program that the second sequence is inverted and assigned to a third sequence. Please, read the first and second sequence from the command line. Calculate the distances between the first and the second and between the first and the third sequence.
Compare the distance between the first and the second and the first and the third sequence and print the alignment with the smaller distance. If the distances are equal, then print the alignment of the first and second sequence.
Test your program with the following sequences:
Solution
import sys
# get sequences from the command line arguments
seq1 = sys.argv[1]
seq2 = sys.argv[2]
### reverse the seq2 string and save as seq2_rev
# initiate variable
seq2_rev = ""
# for each letter in seq2
for i in range(len(seq2)):
# add the next character to the reversed seq2 string
seq2_rev += seq2[len(seq2) - i - 1]
# only run the calculation if the sequences have the same length
if (len(seq1) == len(seq2)):
# initiate counters
dist_1_2 = 0
dist_1_2rev = 0
# for each letter in the sequences
for i in range(len(seq1)):
# if they don't match, add a distance point
if (seq1[i] != seq2[i]):
dist_1_2 += 1
# if they don't match, add a distance point
if (seq1[i] != seq2_rev[i]):
dist_1_2rev += 1
# if the distance seq1seq2 is less or eq to distance seq1seq2_rev
if (dist_1_2 <= dist_1_2rev):
# print the seq1seq2 sequences and distance score
print("Sequence seq1: ", seq1)
print("Sequence seq2: ", seq2)
print("Distance between seq1 and seq2: ", dist_1_2)
# else, if the seq1seq2_rev distance is less than seq1seq2
else:
# print the seq1seq2_rev sequences and distance score
print("Sequence seq1: ", seq1)
print("Sequence seq2_rev: ", seq2_rev)
print("Distance between seq1 and seq2_rev: ", dist_1_2rev)
# tell the user the lengths differ
else:
print("Sequences seq1 and seq2 are of different length.")
Bonus Exercises¶
assignment 2: Functions
-
Open an editor and save your new program. In this program we will create a few functions.
-
Define two functions
similarity
anddistance
:
Note: Purines are A and G, pyrimidines are C and T.
Solution
# define which bases are purines and pyrimidines
pur = ["A", "G"]
pyr = ["C", "T"]
# define the similarity function for two single bases
def similarity(base1, base2):
# if they match, return 1
if (base1 == base2):
return 1
# else,if they dont match but are of the same kind
elif (((base1 in pur) and (base2 in pur)) or ((base1 in pyr) and (base2 in pyr))):
return 0.5
# if they neither matches or are of the same kind, return 0
else:
return 0
# define the distance function for two single bases
def distance(base1, base2):
# if they match, return 0
if (base1 == base2):
return 0
# else,if they dont match but are of the same kind
elif (((base1 in pur) and (base2 in pur)) or ((base1 in pyr) and (base2 in pyr))):
return 0.5
# if they neither matches or are of the same kind, return 1
else:
return 1
- Write two functions
sequence_similarity
andsequence_distance
, which calculates the similarity and distance of two whole sequences.
Solution
# define the similarity function for whole sequences
def sequence_similarity (seq1, seq2):
# initiate counter
similarity_score = 0.0
# go through all bases in seq1
for i in range(len(seq1)):
# calculate their similarity and add to the score
similarity_score = similarity_score + similarity(seq1[i], seq2[i])
# return the final score
return similarity_score
# define the distance function for whole sequences
def sequence_distance(seq1, seq2):
# initiate counter
distance_score = 0.0
# go through all bases in seq1
for i in range(len(seq1)):
# calculate the distance and add to the score
distance_score = distance_score + distance(seq1[i], seq2[i])
# return the final score
return distance_score
- Calculate the similarity and distance for the following sequences. Read these sequences from the command line and print out their similarity and distance.
Solution
import sys
### Paste here the code for the functions you wrote in 1.2 and 1.3 ###
# read the sequences from command line arguments
seq1 = sys.argv[1]
seq2 = sys.argv[2]
# print the similarity and distance
print("Similarity: ", sequence_similarity(seq1, seq2))
print("Distance: ", sequence_distance(seq1, seq2))
assignment 3: Modules
In this exercise we will write three different programs.
-
Write a new Python file (module) called
sequence_tools.py
which contain both the two functionssimilarity
anddistance
as defined previously.Solution
sequence_tools.py######################### ### sequence_tools.py ### ######################### # define which bases are purines and pyrimidines pur = ["A", "G"] pyr = ["C", "T"] # define the similarity function for two single bases def similarity(base1, base2): # if they match, return 1 if (base1 == base2): return 1 # else,if they dont match but are of the same kind elif (((base1 in pur) and (base2 in pur)) or ((base1 in pyr) and (base2 in pyr))) return 0.5 # if they neither matches or are of the same kind, return 0 else: return 0 # define the distance function for two single bases def distance(base1, base2): # if they match, return 0 if (base1 == base2): return 0 # else,if they dont match but are of the same kind elif (((base1 in pur) and (base2 in pur)) or ((base1 in pyr) and (base2 in pyr))) return 0.5 # if they neither matches or are of the same kind, return 1 else: return 1 # define the similarity function for whole sequences def sequence_similarity (seq1, seq2): # initiate counter similarity_score = 0.0 # go through all bases in seq1 for i in range(len(seq1)): # calculate their similarity and add to the score similarity_score = similarity_score + similarity(seq1[i], seq2[i]) # return the final score return similarity_score # define the distance function for whole sequences def sequence_distance(seq1, seq2): # initiate counter distance_score = 0.0 # go through all bases in seq1 for i in range(len(seq1)): # calculate the distance and add to the score distance_score = distance_score + distance(seq1[i], seq2[i]) # return the final score return distance_score
-
Write another Python file that calculates for each combination of two sequences stored in list
seq_list
the similarity and distance using the module defined previously.Solution
main.pyfrom sequence_tools import * # define sequences seq_list = ["ATCCGGT", "GCGTTAC", "CTACTGC", "TTGCAGT", "AGTCACC"] # loop over each sequence in seq_list for i in range(len(seq_list)): # loop over the remaining sequences in seq_list for j in range(i+1, len(seq_list)): # calculate the similarity and distance similarity_score = sequence_similarity(seq_list[i], seq_list[j]) distance_score = sequence_distance(seq_list[i], seq_list[j]) # print the result for this comparison print(seq_list[i], seq_list[j], " Similarity: ", similarity_score, " Distance: ", distance_score)
ATCCGGT GCGTTAC Similarity: 2.5 Distance: 4.5 ATCCGGT CTACTGC Similarity: 3.5 Distance: 3.5 ATCCGGT TTGCAGT Similarity: 4.5 Distance: 2.5 ATCCGGT AGTCACC Similarity: 3.5 Distance: 3.5 GCGTTAC CTACTGC Similarity: 4.0 Distance: 3.0 GCGTTAC TTGCAGT Similarity: 3.0 Distance: 4.0 GCGTTAC AGTCACC Similarity: 2.0 Distance: 5.0 CTACTGC TTGCAGT Similarity: 4.5 Distance: 2.5 CTACTGC AGTCACC Similarity: 2.0 Distance: 5.0 TTGCAGT AGTCACC Similarity: 2.5 Distance: 4.5
-
** Extend your program. Determine the combination of sequences with the highest similarity of all sequences stored in list l. Write these two sequences and the alignment into a new file, called
similar_sequences.txt
.**For example for two given sequences: “ATC” and “ACC” The alignment would be:
And this alignment should be written to a new output file. Hint: A line-break in Python can be made by adding ’\n’ to the end of the line.Solution
main.pyfrom sequence_tools import * # define sequences seq_list = ["ATCCGGT", "GCGTTAC", "CTACTGC", "TTGCAGT", "AGTCACC"] # define variables similarity_highscore = 0 best_seq1 = "" best_seq2 = "" # loop over each sequence in seq_list for i in range(len(seq_list)): # compare the sequence to all remaining sequences in seq_list for j in range(i+1, len(seq_list)): # calculate the similarity similarity_score = sequence_similarity(seq_list[i], seq_list[j]) # check if it's a new similarity highscore if (similarity_score > similarity_highscore): # if it is, save this as the new highscore similarity_highscore = similarity_score best_seq1 = seq_list[i] best_seq2 = seq_list[j] # create an empty string to add the alignment to alignment_matches = "" # go through each letter the best aligned pair for i in range(len(best_seq1)): # find places where they match if (best_seq1[i] == best_seq2[i]): alignment_matches = alignment_matches + "|" # and places they don't match else: alignment_matches = alignment_matches + " " # write the sequences and the match symbols to file outfile = open("similar_sequences.txt", "w") outfile.write(best_seq1 + "\n") outfile.write(alignment_matches + "\n") outfile.write(best_seq2 + "\n")