<jupyter_start><jupyter_text>This is a sandbox/preface for the solubility calculation assignmentThe solubility calculation assignment focuses on building a simple linear solubility model which attempts to predict solubilities for a series of new compounds based on a model trained on a set of compounds with known solubilities. To put it another way, we have a test set and a training set, and want to use the known solubilities from the training set to predict solubilities for the test set.  For solubility prediction, we'll use a series of *descriptors*Descriptors are properties of our molecule which might (or might not) be related to the solubility. For example, we might think that solubility will in general tend to go down as molecular weight goes up, and go up as polarity increases (or go down as polarity decreases) and so on. Here, let's take a sample molecule and calculate a series of descriptors which we might want to use in constructing a simple solubility model.<jupyter_code>from openeye.oechem import *
from openeye.oemolprop import *
from openeye.oeiupac import *
from openeye.oezap import *
from openeye.oeomega import *
import numpy as np
import scipy.stats

#Initialize an OpenEye molecule
mol = OEMol()

#let's look at phenol
OEParseIUPACName( mol, 'naphthalene' )

#Generate conformation
omega = OEOmega()
omega(mol)

#Here one of the descriptors we'll use is the calculated solvation free energy, from OpenEye's ZAP electrostatics solver
#Get zap ready for electrostatics calculations
zap = OEZap()
zap.SetInnerDielectric( 1.0 )
zap.SetGridSpacing(0.5)
area = OEArea()

#Reduce verbosity
OEThrow.SetLevel(OEErrorLevel_Warning)


#Let's print a bunch of properties
#Molecular weight
print( "Molecular weight: %.2f" % OECalculateMolecularWeight(mol) )
#Number of atoms
print( "Number of atoms: %s" % mol.NumAtoms() ) 
#Number of heavy atoms
print( "Number of heavy atoms: %s" % OECount(mol, OEIsHeavy() ) )
#Number of ring atoms
print( "Number of ring atoms: %s" % OECount(mol, OEAtomIsInRing() ) )
#Number of halogens
print( "Number of halogens: %s" % OECount( mol, OEIsHalogen() ))
print ("Number of nitrogens: %s" % OECount( mol, OEIsNitrogen() ) )
print( "Number of oxygens: %s" % OECount( mol, OEIsOxygen() ) )
print( "Number of rotatable bonds: %s" % OECount( mol, OEIsRotor() ) )

#Calculated logP - water to octanol partitioning coefficient (which is often something which may correlate somewhat with solubility)
print( "Calculated logP: %.2f" %  OEGetXLogP( mol ) )

print( "Number of aromatic rings: %s" % OEGetAromaticRingCount( mol ) )

    
    
#Calculate lots of other properties using molprop toolkit as per example in OE MolProp manual
#Handle the setup of 'filter', which computes lots of properties with the goal of filtering compounds. Here we'll not do any filtering
#and will use it solely for property calculation
filt = OEFilter()
ostr = oeosstream()
pwnd = False
filt.SetTable( ostr, pwnd)
#headers = ostr.str().split('\t') #Python 2.x would want something like this; Python 3 version follows
headers = ostr.str().decode().split('\t')
ostr.clear()
filt(mol)
#fields = ostr.str().split('\t') #Python 2.x would want something like this; Python 3 version follows
fields = ostr.str().decode().split('\t')
tmpdct = dict( zip(headers, fields) ) #Format the data we need into a dictionary for easy extraction

print("Polar surface area: %s" % tmpdct[ '2d PSA' ] )
print("Number of hbond donors: %s" % int(tmpdct['hydrogen-bond donors']) )
print("Number of hbond acceptors: %s" % int(tmpdct['hydrogen-bond acceptors']) )
print ("Number of rings: %s" % int(tmpdct['number of ring systems']) )
#print(tmpdct.keys())

#Quickly estimate hydration free energy, or a value correlated with that -- from ZAP manual
#Do ZAP setup for molecule
OEAssignBondiVdWRadii(mol)
OEMMFFAtomTypes(mol)
OEMMFF94PartialCharges(mol)
zap.SetMolecule( mol )
solv = zap.CalcSolvationEnergy()
aval = area.GetArea( mol )
#Empirically estimate solvation free energy (hydration)
solvation = 0.59*solv + 0.01*aval #Convert electrostatic part to kcal/mol; use empirically determined kcal/sq angstrom value times surface area term
print ("Calculated solvation free energy: %.2f" % solvation)<jupyter_output>Molecular weight: 128.17
Number of atoms: 18
Number of heavy atoms: 10
Number of ring atoms: 10
Number of halogens: 0
Number of nitrogens: 0
Number of oxygens: 0
Number of rotatable bonds: 0
Calculated logP: 3.57
Number of aromatic rings: 2
Polar surface area: 0.00
Number of hbond donors: 0
Number of hbond acceptors: 0
Number of rings: 1
Calculated solvation free energy: -4.13<jupyter_text>In the assignment, these get stored in a dictionary. Let's see how that works.<jupyter_code>#Initialize an empty dictionary
compounds = {}
#Name we're working with
molname = 'phenol'
#Create a new OEMol to store this into
mol = OEMol()

#let's look at phenol
OEParseIUPACName( mol, molname )

#Generate conformation
omega = OEOmega()
omega(mol)

#Create a slot in our dictionary for phenol
compounds[molname] = {} #Make it another empty dictionary

#Now let's store some stuff in there
compounds[molname]['mw'] = OECalculateMolecularWeight(mol)
compounds[molname]['rotatable bonds'] = OECount( mol, OEIsRotor() )



#TO DO: Try making an update here to add properties for another compound of your choice to the dictionary



#Let's print it out
print(compounds)<jupyter_output>{'phenol': {'mw': 94.11124000000002, 'rotatable bonds': 0}}<jupyter_text>The point here is just that a dictionary is a flexible data structure which allows us to easily store away information we might want later in an organized manner. For example, if I want to see everything I have for phenol, I simply use:<jupyter_code>print( compounds['phenol'])<jupyter_output>{'mw': 94.11124000000002, 'rotatable bonds': 0}
{'fraction_non_alphanumeric': 0.08931034482758621, 'fraction_numerical': 0.014137931034482758, 'mean_word_length': 5.284940411700975, 'pattern_counts': {'":': 0, '<': 10, '<?xml version=': 0, '>': 10, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 2, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
{'chain_length': '3', 'dir': 'jupyter-structured-clean-dedup', 'id': '105659', 'license': 'CC-BY-4.0', 'n_tokens_mistral': 1955, 'n_tokens_neox': 1739, 'n_words': 749, 'path': 'uci-pharmsci/lectures/empirical_physical_properties/physprops_solubility.ipynb', 'repo_name': 'davidlmobley/drug-computing', 'size': '5800'}