
   0 | ||| This computes the logP in accordance with the JLogP module in
   1 | ||| the CDK.
   2 | |||
   3 | ||| This was published in the
   4 | ||| [J. Chem. Inf.](https://jcheminf.biomedcentral.com/articles/10.1186/s13321-018-0316-5)
   5 | module Chem.QSAR.JPlogP
   6 | 
   7 | import Chem
   8 | import Chem.Aromaticity
   9 | import Chem.QSAR.Util
  10 | import Data.SortedMap
  11 | 
  12 | %default total
  13 | 
  14 | --------------------------------------------------------------------------------
  15 | -- Utilities
  16 | --------------------------------------------------------------------------------
  17 | 
  18 | isPolar : Elem -> Bool
  19 | isPolar O = True
  20 | isPolar N = True
  21 | isPolar S = True
  22 | isPolar P = True
  23 | isPolar _ = False
  24 | 
  25 | isElectronWithdrawing : Elem -> Bool
  26 | isElectronWithdrawing O  = True
  27 | isElectronWithdrawing N  = True
  28 | isElectronWithdrawing S  = True
  29 | isElectronWithdrawing F  = True
  30 | isElectronWithdrawing Cl = True
  31 | isElectronWithdrawing Br = True
  32 | isElectronWithdrawing I  = True
  33 | isElectronWithdrawing _  = False
  34 | 
  35 | -- fast failing addition of floating point values with the potential of failure
  36 | sumMaybes : Double -> (a -> Maybe Double) -> List a -> Maybe Double
  37 | sumMaybes v f []      = Just v
  38 | sumMaybes v f (x::xs) =
  39 |   case f x of
  40 |     Just y  => sumMaybes (v+y) f xs
  41 |     Nothing => Nothing
  42 | 
  43 | -- list of pairs at the bottom of the source file
  44 | weights : SortedMap Nat Double
  45 | 
  46 | parameters {0 b,e,p,r,t,c,l : Type}
  47 |            {auto ce : Cast e Elem}
  48 |            {auto cb : Cast b BondOrder}
  49 |            {k       : Nat}
  50 |            (g       : IGraph k (AromBond b) (Atom e Charge p r HCount t c l))
  51 | 
  52 |   subsel : Elem -> HCount -> AssocList k (AromBond b) -> Nat
  53 | 
  54 |   hydrogenSpecial : List (Fin k) -> Nat
  55 | 
  56 |   hydrogenCode : Fin k -> Nat
  57 |   hydrogenCode n = 101_100 + hydrogenSpecial [n]
  58 | 
  59 |   -- returns an atom code in the form CAANSS where
  60 |   -- C  = charge+1
  61 |   -- AA = Atomic Number
  62 |   -- N  = NonH Neighbour count
  63 |   -- SS = elemental subselection via bonding and neighbour identification
  64 |   atomCode : Fin k -> Nat
  65 |   atomCode n =
  66 |     let A a ns := adj g n
  67 |         el     := cast @{ce} (elem a)
  68 |         chrg   := a.charge.value + 1
  69 |         nonHs  := countNonHs g ns
  70 |         anum   := cast {to = Nat} $ value $ atomicNr el
  71 |         add    := subsel el a.hydrogen ns
  72 |      in if add > 99 || anum > 99 || nonHs > 9 || chrg < 0
  73 |            then 0
  74 |            else cast chrg * 100_000 + anum * 1000 + nonHs * 100 + add
  75 | 
  76 |   ||| Computes the logP contribution for a single atom and its
  77 |   ||| implicit hydrogen atoms.
  78 |   |||
  79 |   ||| Returns `Nothing` in case the atom is of an unsupported element
  80 |   ||| or connectivity.
  81 |   export
  82 |   jpLogPContrib : Fin k -> Maybe Double
  83 |   jpLogPContrib n =
  84 |     case hydrogen $ lab g n of
  85 |       0 => lookup (atomCode n) weights
  86 |       h =>
  87 |         let Just cn := lookup (atomCode n) weights     | _ => Nothing
  88 |             Just ch := lookup (hydrogenCode n) weights | _ => Nothing
  89 |          in Just $ cast h.value * ch + cn
  90 |   
  91 |   ||| Sums up the logP contributions for all atoms in a molecule.
  92 |   |||
  93 |   ||| Returns `Nothing` in case `jpLogPContrib` returns `Nothing` for
  94 |   ||| at least one atom.
  95 |   export %inline
  96 |   jpLogP : Maybe Double
  97 |   jpLogP = sumMaybes 0.0 jpLogPContrib (nodes g)
  98 | 
  99 | --------------------------------------------------------------------------------
 100 | -- Implementation
 101 | --------------------------------------------------------------------------------
 102 | 
 103 |   addBond : (Elem -> Bool) -> AromBonds -> (Fin k, AromBond b) -> AromBonds
 104 |   addBond p bs (n,AB t ar) =
 105 |     let True  := p (elemAt g n) | False => bs
 106 |         False := ar             | True  => {arom $= S} bs
 107 |      in case cast @{cb} t of
 108 |           Single => {single $= S} bs
 109 |           Dbl    => {double $= S} bs
 110 |           Triple => {triple $= S} bs
 111 | 
 112 |   bondsWhere : (Elem -> Bool) -> AssocList k (AromBond b) -> AromBonds
 113 |   bondsWhere p = foldl (addBond p) (ABS 0 0 0 0) . pairs
 114 | 
 115 |   bondCount : Fin k -> Nat
 116 |   bondCount n = let A a ns := adj g n in numNeighbours a ns
 117 | 
 118 |   -- formal oxidation state computed from the bond order of
 119 |   -- electron withdrawing neighbours.
 120 |   formalOx : Fin k -> Nat
 121 |   formalOx = foldl acc 0 . neighboursAsPairs g
 122 |     where
 123 |       acc : Nat -> (Fin k, AromBond b) -> Nat
 124 |       acc c (n,b) =
 125 |         case (isElectronWithdrawing (elemAt g n), cast @{cb} b.type) of
 126 |           (True, Single) => c + 1
 127 |           (True, Dbl)    => c + 2
 128 |           (True, Triple) => c + 3
 129 |           _              => c
 130 | 
 131 |   alphaCarbonyl : Elem -> AssocList k (AromBond b) -> Bool
 132 |   alphaCarbonyl el = any (hasNeighbour g el Single) . keys
 133 | 
 134 |   nextToAromatic : AssocList k (AromBond b) -> Bool
 135 |   nextToAromatic ns =
 136 |     not (any AromBond.arom ns) &&
 137 |     any (\(n,b) => Single == cast b && isAromatic g n) (pairs ns)
 138 | 
 139 |   doubleBondHetero : AssocList k (AromBond b) -> Bool
 140 |   doubleBondHetero =
 141 |     any (\(n,AB t b) => not b && cast t == Dbl && isPolar (elemAt g n)) . pairs
 142 | 
 143 |   isCarbonyl : (Fin k, AromBond b) -> Bool
 144 |   isCarbonyl (n, AB t b) =
 145 |     not b && cast t == Single && doubleBondHetero (neighboursAsAL g n)
 146 | 
 147 |   fluorineSpecial : List (Fin k) -> Nat
 148 |   fluorineSpecial [n] =
 149 |     case elemAt g n of
 150 |       S => 8 -- F-S
 151 |       B => 9 -- F-B
 152 |       C => case bondCount n of
 153 |         2 => 2
 154 |         3 => 3
 155 |         4 => if formalOx n <= 2 then 5 else 7
 156 |         _ => 99
 157 |       _ => 1 -- F-hetero
 158 |   fluorineSpecial _  = 99
 159 | 
 160 |   oxygenSpecial : HCount -> AssocList k (AromBond b) -> Nat
 161 |   oxygenSpecial h ns =
 162 |     case cast h.value + length ns of
 163 |       2 =>
 164 |         if      hasElem g N ns then 1
 165 |         else if hasElem g S ns then 2
 166 |         else if any AromBond.arom ns then 8
 167 |         else    3
 168 |       1 =>
 169 |         if      hasElem g N ns then 4
 170 |         else if hasElem g S ns then 5
 171 |         else if alphaCarbonyl O ns then 6
 172 |         else if alphaCarbonyl N ns then 9
 173 |         else if alphaCarbonyl S ns then 10
 174 |         else    7
 175 |       _ => 0
 176 | 
 177 |   nitrogenSpecial : HCount -> AssocList k (AromBond b) -> Nat
 178 |   nitrogenSpecial h ns =
 179 |     case cast h.value + length ns of
 180 |       4 => 9
 181 |       3 =>
 182 |         if      nextToAromatic ns                  then 1
 183 |         else if any isCarbonyl (pairs ns)          then 2
 184 |         else if doubleBondHetero ns                then 10
 185 |         else if single (bondsWhere isPolar ns) > Z then 3
 186 |         else                                 4
 187 |       2 =>
 188 |         if      any AromBond.arom ns then 5
 189 |         else if doubleBondHetero ns  then 6
 190 |         else                              7
 191 |       1 => 8
 192 |       _ => 0
 193 | 
 194 |   carbonSpecial : HCount -> AssocList k (AromBond b) -> Nat
 195 |   carbonSpecial h ns =
 196 |     case cast h.value + length ns of
 197 |       4 => if single (bondsWhere isPolar ns) > Z then 3 else 2
 198 |       3 => case any AromBond.arom ns of
 199 |         True => case bondsWhere isPolar ns of
 200 |           ABS 0 (S _) _ _ => 11
 201 |           ABS 1 0     _ _ => 5
 202 |           ABS 1 _     _ _ => 13
 203 |           _              => 4
 204 |         False => case bondsWhere isPolar ns of
 205 |           ABS 0     _ 1 _ => 7
 206 |           ABS _     _ 1 _ => 14
 207 |           ABS (S _) _ 0 _ => 8
 208 |           _              => 6
 209 |       2 => case bondsWhere isPolar ns of
 210 |         ABS 0 _ _ 1 => 12
 211 |         ABS 1 _ _ 0 => 10
 212 |         ABS 1 _ _ 1 => 15
 213 |         _          => 9
 214 |       _ => if numBonds (bondsWhere isPolar ns) > Z then 1 else 0
 215 | 
 216 |   hydrogenSpecial [n] =
 217 |     let A a ns := adj g n
 218 |      in case cast @{ce} (elem a) of
 219 |           C => case any isCarbonyl (pairs ns) of
 220 |             True  => 51
 221 |             False => case numNeighbours a ns of
 222 |               4 => case formalOx n of
 223 |                 0 => 46
 224 |                 1 => 47
 225 |                 2 => 48
 226 |                 _ => 49
 227 |               3 => case formalOx n of
 228 |                 0 => 47
 229 |                 1 => 48
 230 |                 _ => 49
 231 |               2 => case formalOx n of
 232 |                 0 => 48
 233 |                 _ => 49
 234 |               1 => 121
 235 |               _ => 0
 236 |           _ => 50
 237 |   hydrogenSpecial _   = 0
 238 | 
 239 |   subsel C h ns = carbonSpecial h ns
 240 |   subsel N h ns = nitrogenSpecial h ns
 241 |   subsel O h ns = oxygenSpecial h ns
 242 |   subsel H h ns = hydrogenSpecial (keys ns)
 243 |   subsel F h ns = fluorineSpecial (keys ns)
 244 |   subsel _ h ns =
 245 |     if any AromBond.arom ns then 10 else numBonds $ bondsWhere isPolar ns
 246 | -- -- 
 247 | -- -- 
 248 | 
 249 | weights =
 250 |   fromList
 251 |     [ (115201, 0.09994758075256505)
 252 |     , (115200, 0.3379378715836258)
 253 |     , (134400, -0.601185704550091)
 254 |     , (115202, 0.30788026393512663)
 255 |     , (207110, -0.26496784659823264)
 256 |     , (134404, -1.1724223083800398)
 257 |     , (115210, -0.08346526510422402)
 258 |     , (153100, 0.9270429695306335)
 259 |     , (153101, 1.0145354986151272)
 260 |     , (116503, 0.4425591506257104)
 261 |     , (133402, -0.2557512835716269)
 262 |     , (114200, 0.01526633068977459)
 263 |     , (133403, -0.8169297847709985)
 264 |     , (5401, 0.10441747048024147)
 265 |     , (101147, 1.2616122128062803)
 266 |     , (5402, 0.05162677089603265)
 267 |     , (101146, 1.3994445700193028)
 268 |     , (133401, -0.3639701318790265)
 269 |     , (5403, 0.6788714142147848)
 270 |     , (101149, 1.3258747052968567)
 271 |     , (101148, 1.2711079053599976)
 272 |     , (101151, 1.2556350911435799)
 273 |     , (133404, -0.6891007636859257)
 274 |     , (101150, -0.2363827296335956)
 275 |     , (107301, 0.1787473725640726)
 276 |     , (107303, -0.016959741231455404)
 277 |     , (107302, 0.17510323694483412)
 278 |     , (7207, -0.1902718970453204)
 279 |     , (107304, 0.1542614658530209)
 280 |     , (109101, 0.41374364339440817)
 281 |     , (115501, -0.14068629542864905)
 282 |     , (115500, 0.17750686328369028)
 283 |     , (115503, 0.013887172778612027)
 284 |     , (109103, 0.26651823980406203)
 285 |     , (115502, -0.11992335751384754)
 286 |     , (115505, -0.34311166776744884)
 287 |     , (109105, 0.43019405241170144)
 288 |     , (115504, -0.1025811855926768)
 289 |     , (207409, -0.23852255872700964)
 290 |     , (109107, 0.46487210540519147)
 291 |     , (109109, 0.4801727828186138)
 292 |     , (109108, 0.2430918227919903)
 293 |     , (134202, -0.631693669606887)
 294 |     , (134200, -0.04910253697266963)
 295 |     , (134201, 0.011171177597208612)
 296 |     , (106303, -1.5239332211866237)
 297 |     , (106302, -1.3023723757750838)
 298 |     , (106305, 0.11154050989797104)
 299 |     , (134210, 0.8725452294362313)
 300 |     , (106304, 0.20930194677993302)
 301 |     , (106307, -0.2690238196488019)
 302 |     , (106306, 0.01713355115342431)
 303 |     , (108101, -0.09994092418236927)
 304 |     , (106308, 0.026068238656834854)
 305 |     , (108103, 0.009836423634389751)
 306 |     , (106311, 0.1427758223243856)
 307 |     , (108102, 0.09143879448168145)
 308 |     , (108105, -0.4749180976677123)
 309 |     , (106313, 0.2897701159093737)
 310 |     , (108104, 0.015474674863268114)
 311 |     , (108107, 0.03602987937616161)
 312 |     , (108106, 0.19048034389371205)
 313 |     , (106314, -0.39178934954486583)
 314 |     , (108109, -0.10666832936758427)
 315 |     , (116301, -0.04914876490157527)
 316 |     , (116300, 0.7367961788572792)
 317 |     , (116303, -0.16169601215900375)
 318 |     , (116302, -0.12643665926860584)
 319 |     , (108110, 0.08928780909185945)
 320 |     , (5201, 0.2201279736215536)
 321 |     , (5202, 0.19045980382858146)
 322 |     , (133200, -0.3076946592117401)
 323 |     , (208208, 0.7099234399396344)
 324 |     , (133201, -0.49680932374826)
 325 |     , (105301, -0.1621399782797916)
 326 |     , (105300, -0.174370011345452)
 327 |     , (105303, -0.1432571497021001)
 328 |     , (105302, -0.17755682989875274)
 329 |     , (107101, 0.051702151541644426)
 330 |     , (107103, -0.2691841786263712)
 331 |     , (107102, 0.06496457627779738)
 332 |     , (107104, -0.33802382979998147)
 333 |     , (107107, 0.394978253857827)
 334 |     , (107106, 0.3859974866654674)
 335 |     , (207207, -0.07239795705976523)
 336 |     , (115301, -0.023836916386805847)
 337 |     , (107108, 0.11642255395641928)
 338 |     , (207206, -0.24558051744952064)
 339 |     , (115300, 0.08797456644925557)
 340 |     , (115303, -0.23605983536956895)
 341 |     , (115302, -0.10814292962539623)
 342 |     , (117101, 0.7369857420763359)
 343 |     , (117100, 0.7116079866622599)
 344 |     , (153200, 0.7888787630537003)
 345 |     , (106103, -3.767115237033892)
 346 |     , (106102, -3.616478490407742)
 347 |     , (116600, 1.2424471654297324)
 348 |     , (106107, -2.416958126564593)
 349 |     , (106106, -2.0565095206356196)
 350 |     , (114301, -0.13761744158191205)
 351 |     , (106109, -0.929959267287108)
 352 |     , (114300, -0.3058393642193975)
 353 |     , (114302, -0.3095457315739295)
 354 |     , (106112, -1.3012751020335893)
 355 |     , (116101, 1.2551746199963494)
 356 |     , (116100, 0.7001698255404422)
 357 |     , (105100, 0.36881886842575007)
 358 |     , (216210, 0.7113783097640652)
 359 |     , (134302, -0.37554769340770927)
 360 |     , (134303, -0.36036185997589393)
 361 |     , (115100, 0.6096777283013224)
 362 |     , (134300, -0.4657894122488925)
 363 |     , (134301, -0.3795150596315356)
 364 |     , (106403, -0.6035455702256183)
 365 |     , (106402, -0.19123067665076543)
 366 |     , (8104, -0.24195926611016633)
 367 |     , (108201, 0.030702954811754002)
 368 |     , (8105, -0.4215369017701643)
 369 |     , (108203, 0.16547595574733062)
 370 |     , (8106, -0.2964579842443157)
 371 |     , (108202, 0.12058552604519218)
 372 |     , (116401, 0.7460081218102875)
 373 |     , (116400, 0.9078902724309305)
 374 |     , (108208, 0.06665724849079398)
 375 |     , (116403, 0.3132146936243478)
 376 |     , (116402, 0.5536499705080733)
 377 |     , (133302, -0.7030388263360682)
 378 |     , (114100, 0.3587262776601395)
 379 |     , (133303, -0.4317778955330324)
 380 |     , (116404, 0.0372737885433136)
 381 |     , (133300, -0.2042709645502799)
 382 |     , (133301, 0.25473886515880656)
 383 |     , (135100, 0.8603892414982898)
 384 |     , (135101, 0.7235643505309818)
 385 |     , (107201, 0.25454958160787483)
 386 |     , (107203, 0.20036994532453556)
 387 |     , (107202, 0.15973823557158393)
 388 |     , (107205, -0.2543417785988061)
 389 |     , (7108, 0.2634019114789942)
 390 |     , (207303, -0.47290811359153156)
 391 |     , (107204, 0.1275311923081761)
 392 |     , (207302, 0.24535814336356004)
 393 |     , (107207, 0.3145047278208787)
 394 |     , (207301, 0.36070459894156437)
 395 |     , (107206, 0.3113299058370724)
 396 |     , (115401, -0.33332080419904164)
 397 |     , (115400, -0.10861521692106911)
 398 |     , (115403, -0.7356011823089021)
 399 |     , (207304, -0.07464529856367844)
 400 |     , (115402, -0.37912343134016635)
 401 |     , (115404, -0.8103937160471404)
 402 |     , (207310, 0.35707776316923207)
 403 |     , (153302, 0.6015475892374368)
 404 |     , (134100, -0.1363315394071071)
 405 |     , (134101, 0.1092440266426836)
 406 |     , (153301, 0.5748646504612611)
 407 |     , (106203, -2.6563042881537227)
 408 |     , (106202, -2.349643420766774)
 409 |     , (106204, -0.953850490807928)
 410 |     , (106207, -1.5339493033443787)
 411 |     , (106206, -1.0109064714897187)
 412 |     , (106209, 0.19490613802180773)
 413 |     , (114401, -0.16389842380630032)
 414 |     , (106208, -1.195325147036)
 415 |     , (114400, -0.11455385699777243)
 416 |     , (106211, -1.1770219185049515)
 417 |     , (114403, -0.1553108949228919)
 418 |     , (114402, -0.1103226433176957)
 419 |     , (106210, 0.16163164657523407)
 420 |     , (106212, -0.2512600932013654)
 421 |     , (114404, -0.13111265719465162)
 422 |     , (106215, -0.14610561930731214)
 423 |     , (106214, -1.6084019086078098)
 424 |     , (116201, 0.6390448917281739)
 425 |     , (116200, 0.7928741182178681)
 426 |     , (116202, 0.5335850658686033)
 427 |     , (216301, 0.16768140357444056)
 428 |     , (216300, 0.8104414861979382)
 429 |     , (105201, 0.06904337130830694)
 430 |     , (105202, -0.0745914491562598)
 431 |     , (116210, 0.5791544003852439)
 432 |     , (216310, 0.16982252294512776)
 433 |     ]
 434 |