diff --git a/propka/atom.py b/propka/atom.py
index 3fbaa86..40e15b2 100644
--- a/propka/atom.py
+++ b/propka/atom.py
@@ -104,7 +104,7 @@ class Atom(object):
             self.y = float(line[38:46].strip())
             self.z = float(line[46:54].strip())
             self.res_num = int(line[22:26].strip())
-            self.res_name = "%-3s" % (line[17:20].strip())
+            self.res_name = "{0:<3s}".format(line[17:20].strip())
             self.chain_id = line[21]
             # Set chain id to "_" if it is just white space.
             if not self.chain_id.strip():
@@ -124,7 +124,7 @@ class Atom(object):
             if len(self.name) == 4:
                 self.element = self.element[0]
             if len(self.element) == 2:
-                self.element = '%1s%1s' % (
+                self.element = '{0:1s}{1:1s}'.format(
                     self.element[0], self.element[1].lower())
 
     def set_group_type(self, type_):
@@ -278,15 +278,15 @@ class Atom(object):
         Returns:
             String with PDB line.
         """
-        res = 'CONECT%5d' % self.numb
+        res = 'CONECT{0:5d}'.format(self.numb)
 
         bonded = []
         for atom in self.bonded_atoms:
             bonded.append(atom.numb)
         bonded.sort()
 
-        for b in bonded:
-            res += '%5d' % b
+        for bond in bonded:
+            res += '{0:5d}'.format(bond)
         res += '\n'
         return res
 
@@ -316,12 +316,14 @@ class Atom(object):
             self.occ = self.occ.replace('LG', 'non_titratable_ligand')
             # try to initialise the group
             try:
-                group_attr = "%s_group" % self.occ
+                group_attr = "{0:s}_group".format(self.occ)
                 group_attr = getattr(propka.group, group_attr)
                 self.group = group_attr(self)
             except:
                 # TODO - be more specific with expection handling here
-                str_ = '%s in input_file is not recognized as a group' % self.occ
+                str_ = (
+                    '{0:s} in input_file is not recognized as a group'.format(
+                        self.occ))
                 raise Exception(str_)
         # set the model pKa value
         if self.beta != '-':
diff --git a/propka/bonds.py b/propka/bonds.py
index 0144f78..c6656df 100644
--- a/propka/bonds.py
+++ b/propka/bonds.py
@@ -191,7 +191,7 @@ class BondMaker:
             atoms:  list of atoms to check for bonds
         """
         for atom1 in atoms:
-            key = '%s-%s' % (atom1.res_name, atom1.name)
+            key = '{0:s}-{1:s}'.format(atom1.res_name, atom1.name)
             if key in list(self.num_pi_elec_bonds_sidechains.keys()):
                 atom1.num_pi_elec_2_3_bonds = (
                     self.num_pi_elec_bonds_sidechains[key])
@@ -235,7 +235,7 @@ class BondMaker:
                         self.num_pi_elec_conj_bonds_ligands[atom.sybyl_type])
             # for protein
             if atom.type == 'atom':
-                key = '%s-%s' % (atom.res_name, atom.name)
+                key = '{0:s}-{1:s}'.format(atom.res_name, atom.name)
                 if key in list(self.num_pi_elec_bonds_sidechains.keys()):
                     atom.num_pi_elec_2_3_bonds = (
                         self.num_pi_elec_bonds_sidechains[key])
@@ -300,7 +300,7 @@ class BondMaker:
         sq_dist = propka.calculations.squared_distance(atom1, atom2)
         if sq_dist > self.max_sq_distance:
             return False
-        key = '%s-%s' % (atom1.element, atom2.element)
+        key = '{0:s}-{1:s}'.format(atom1.element, atom2.element)
         h_count = key.count('H')
         if sq_dist < self.h_dist_squared and h_count == 1:
             return True
diff --git a/propka/calculations.py b/propka/calculations.py
index c7c84af..3770983 100644
--- a/propka/calculations.py
+++ b/propka/calculations.py
@@ -242,9 +242,9 @@ def add_trp_hydrogen(residue):
         elif atom.name == "CE2":
             ce_atom = atom
     if (cd_atom is None) or (ne_atom is None) or (ce_atom is None):
-        errstr = "Unable to find all atoms for %s %s" % (
+        str_ = "Unable to find all atoms for {0:s} {1:s}".format(
             residue[0].res_name, residue[0].res_num)
-        raise ValueError(errstr)
+        raise ValueError(str_)
     he_atom = protonate_sp2(cd_atom, ne_atom, ce_atom)
     he_atom.name = "HNE"
 
@@ -269,9 +269,9 @@ def add_amd_hydrogen(residue):
               or (atom.res_name == "ASN" and atom.name == "ND2")):
             n_atom = atom
     if (c_atom is None) or (o_atom is None) or (n_atom is None):
-        errstr = "Unable to find all atoms for %s %s" % (
+        str_ = "Unable to find all atoms for {0:s} {1:s}".format(
             residue[0].res_name, residue[0].res_num)
-        raise ValueError(errstr)
+        raise ValueError(str_)
     h1_atom = protonate_direction(n_atom, o_atom, c_atom)
     h1_atom.name = "HN1"
     h2_atom = protonate_average_direction(n_atom, c_atom, o_atom)
@@ -396,10 +396,10 @@ def make_new_h(atom, x, y, z):
         new hydrogen atom
     """
     new_h = propka.atom.Atom()
-    new_h.set_property(numb=None, name='H%s' % atom.name[1:],
-                       res_name=atom.res_name, chain_id=atom.chain_id,
-                       res_num=atom.res_num, x=x, y=y, z=z, occ=None,
-                       beta=None)
+    new_h.set_property(
+        numb=None, name='H{0:s}'.format(atom.name[1:]),
+        res_name=atom.res_name, chain_id=atom.chain_id,
+        res_num=atom.res_num, x=x, y=y, z=z, occ=None, beta=None)
     new_h.element = 'H'
     new_h.bonded_atoms = [atom]
     new_h.charge = 0
@@ -604,8 +604,9 @@ def hydrogen_bond_interaction(group1, group2, version):
     atoms2 = group2.get_interaction_atoms(group1)
     [closest_atom1, dist, closest_atom2] = get_smallest_distance(atoms1, atoms2)
     if None in [closest_atom1, closest_atom2]:
-        warning('Side chain interaction failed for %s and %s' % (
-            group1.label, group2.label))
+        warning(
+            'Side chain interaction failed for {0:s} and {1:s}'.format(
+                group1.label, group2.label))
         return None
     # get the parameters
     [dpka_max, cutoff] = version.get_hydrogen_bond_parameters(closest_atom1,
diff --git a/propka/conformation_container.py b/propka/conformation_container.py
index f86b739..d44e30c 100644
--- a/propka/conformation_container.py
+++ b/propka/conformation_container.py
@@ -58,7 +58,7 @@ class ConformationContainer:
         # if a group is coupled and we are reading a .propka_input file, then
         # some more configuration might be needed
         map_ = make_interaction_map(
-            'Covalent coupling map for %s' % self,
+            'Covalent coupling map for {0:s}'.format(str(self)),
             self.get_covalently_coupled_groups(),
             lambda g1, g2: g1 in g2.covalently_coupled_groups)
         info(map_)
@@ -100,7 +100,7 @@ class ConformationContainer:
             self.set_common_charge_centres()
         # print coupling map
         map_ = make_interaction_map(
-            'Covalent coupling map for %s' % self,
+            'Covalent coupling map for {0:s}'.format(str(self)),
             self.get_covalently_coupled_groups(),
             lambda g1, g2: g1 in g2.covalently_coupled_groups)
         info(map_)
diff --git a/propka/coupled_groups.py b/propka/coupled_groups.py
index 942fed8..aeed93e 100644
--- a/propka/coupled_groups.py
+++ b/propka/coupled_groups.py
@@ -193,7 +193,7 @@ class NonCovalentlyCoupledGroups:
             conformation:  conformation to print
         """
         map_ = make_interaction_map(
-            'Non-covalent coupling map for %s' % conformation,
+            'Non-covalent coupling map for {0:s}'.format(str(conformation)),
             conformation.get_non_covalently_coupled_groups(),
             lambda g1, g2: g1 in g2.non_covalently_coupled_groups)
         info(map_)
@@ -209,7 +209,8 @@ class NonCovalentlyCoupledGroups:
             conformation:  conformation to print
             system:  system to print
         """
-        info('System containing %d groups:' % len(system))
+        info(
+            'System containing {0:d} groups:'.format(len(system)))
         # make list of interactions within this system
         interactions = list(itertools.combinations(system, 2))
         # print out coupling info for each interaction
@@ -233,8 +234,8 @@ class NonCovalentlyCoupledGroups:
             # Tell the user what is swap in this combination
             swap_info += 'Swapping the following interactions:\n'
             for interaction in combination:
-                swap_info += ' %s %s\n' % (interaction[0].label,
-                                           interaction[1].label)
+                swap_info += ' {0:s} {1:s}\n'.format(
+                    interaction[0].label, interaction[1].label)
             # swap...
             for interaction in combination:
                 self.swap_interactions([interaction[0]], [interaction[1]])
@@ -275,9 +276,9 @@ class NonCovalentlyCoupledGroups:
         all_labels = [g.label for g in system]
         str_ = ' ' + '-' * 113 + '\n'
         for group in system:
-            str_ += self.tagged_format(' %-8s|' % tag,
-                                       group.get_determinant_string(),
-                                       all_labels)
+            str_ += self.tagged_format(
+                ' {0:<8s}|'.format(tag), group.get_determinant_string(),
+                all_labels)
         return str_ + '\n'
 
     def swap_interactions(self, groups1, groups2, include_side_chain_hbs=True):
@@ -342,10 +343,10 @@ class NonCovalentlyCoupledGroups:
         Returns:
             tagged string
         """
-        str_ = "%s %s" % (tag, str_)
-        str_ = str_.replace('\n', '\n%s ' % tag)
+        str_ = "{0:s} {1:s}".format(tag, str_)
+        str_ = str_.replace('\n', '\n{0:s} '.format(tag))
         for label in labels:
-            str_ = str_.replace(label, '\033[31m%s\033[30m' % label)
+            str_ = str_.replace(label, '\033[31m{0:s}\033[30m'.format(label))
         return str_ + '\n'
 
     @staticmethod
diff --git a/propka/determinant.py b/propka/determinant.py
index 2b29d08..5370798 100644
--- a/propka/determinant.py
+++ b/propka/determinant.py
@@ -34,4 +34,4 @@ class Determinant:
         self.value += value
 
     def __str__(self):
-        return '%s: %8.2f' % (self.label, self.value)
+        return '{0:s}: {1:8.2f}'.format(self.label, self.value)
diff --git a/propka/group.py b/propka/group.py
index 9b5c445..125b0dc 100644
--- a/propka/group.py
+++ b/propka/group.py
@@ -179,7 +179,7 @@ class Group:
         # first check if there are any coupled groups at all
         if len(self.covalently_coupled_groups) == 0:
             return ''
-        line = 'CCOUPL%5d' % self.atom.numb
+        line = 'CCOUPL{0:5d}'.format(self.atom.numb)
         # extract and sort numbers of coupled groups
         coupled = []
         for group in self.covalently_coupled_groups:
@@ -187,7 +187,7 @@ class Group:
         coupled.sort()
         # write 'em out
         for num in coupled:
-            line += '%5d' % num
+            line += '{0:5d}'.format(num)
         line += '\n'
         return line
 
@@ -200,7 +200,7 @@ class Group:
         # first check if there are any coupled groups at all
         if len(self.non_covalently_coupled_groups) == 0:
             return ''
-        line = 'NCOUPL%5d' % self.atom.numb
+        line = 'NCOUPL{0:5d}'.format(self.atom.numb)
         # extract and sort numbers of coupled groups
         coupled = []
         for group in self.non_covalently_coupled_groups:
@@ -208,7 +208,7 @@ class Group:
         coupled.sort()
         # write 'em out
         for num in coupled:
-            line += '%5d' % num
+            line += '{0:5d}'.format(num)
         line += '\n'
         return line
 
@@ -229,9 +229,10 @@ class Group:
 
     def __iadd__(self, other):
         if self.type != other.type:
-            errstr = ('Cannot add groups of different types (%s and %s)'
-                      % (self.type, other.type))
-            raise Exception(errstr)
+            str_ = (
+                'Cannot add groups of different types '
+                '({0:s} and {1:s})'.format(self.type, other.type))
+            raise Exception(str_)
         # add all values
         self.pka_value += other.pka_value
         self.num_volume += other.num_volume
@@ -343,9 +344,9 @@ class Group:
             if not self.model_pka_set:
                 self.model_pka = self.parameters.model_pkas[self.residue_type]
                 # check if we should apply a custom model pka
-                key = (
-                    '%s-%s'
-                    % (self.atom.res_name.strip(), self.atom.name.strip()))
+                key = '{0:s}-{1:s}'.format(
+                    self.atom.res_name.strip(),
+                    self.atom.name.strip())
                 if key in self.parameters.custom_model_pkas.keys():
                     self.model_pka = self.parameters.custom_model_pkas[key]
                 self.model_pka_set = True
@@ -388,10 +389,10 @@ class Group:
                         ok = False
         if not ok:
             str_ = 'Missing atoms or failed protonation for '
-            str_ += ('%s (%s) -- please check the structure'
-                     % (self.label, self.type))
+            str_ += '{0:s} ({1:s}) -- please check the structure'.format(
+                self.label, self.type)
             warning(str_)
-            warning('%s' % self)
+            warning('{0:s}'.format(str(self)))
             num_acid = sum(
                 [EXPECTED_ATOMS_ACID_INTERACTIONS[self.type][e]
                  for e in EXPECTED_ATOMS_ACID_INTERACTIONS[self.type].keys()])
@@ -399,15 +400,19 @@ class Group:
                 [EXPECTED_ATOMS_BASE_INTERACTIONS[self.type][e]
                  for e in EXPECTED_ATOMS_BASE_INTERACTIONS[self.type].keys()])
             warning(
-                'Expected %d interaction atoms for acids, found:' % num_acid)
+                'Expected {0:d} interaction atoms for acids, found:'.format(
+                    num_acid))
             for i in range(len(self.interaction_atoms_for_acids)):
                 warning(
-                    '             %s' % self.interaction_atoms_for_acids[i])
+                    '             {0:s}'.format(
+                        str(self.interaction_atoms_for_acids[i])))
             warning(
-                'Expected %d interaction atoms for bases, found:' % num_base)
+                'Expected {0:d} interaction atoms for bases, found:'.format(
+                    num_base))
             for i in range(len(self.interaction_atoms_for_bases)):
                 warning(
-                    '             %s' % self.interaction_atoms_for_bases[i])
+                    '             {0:s}'.format(
+                        str(self.interaction_atoms_for_bases[i])))
 
     def get_interaction_atoms(self, interacting_group):
         """Get atoms involved in interaction with other group.
@@ -461,14 +466,14 @@ class Group:
                               number_of_coulomb)
         str_ = ""
         for line_number in range(number_of_lines):
-            str_ += "%s" % (self.label)
+            str_ += "{0:s}".format(self.label)
             if line_number == 0:
-                str_ += " %6.2lf" %(self.pka_value)
+                str_ += " {0:6.2f}".format(self.pka_value)
                 if len(self.non_covalently_coupled_groups) > 0:
                     str_ += '*'
                 else:
                     str_ += ' '
-                str_ += " %4d%2s " % (int(100.0*self.buried), "%")
+                str_ += " {0:4d}{1:>2s} ".format(int(100.0*self.buried), "%")
                 str_ += " %6.2lf %4d" % (self.energy_volume, self.num_volume)
                 str_ += " %6.2lf %4d" % (self.energy_local, self.num_local)
             else: