0 | module Chem.QSAR.TPSA
  1 |
  2 | import Chem
  3 | import Chem.Aromaticity
  4 | import Chem.AtomType
  5 | import Chem.QSAR.Util
  6 | import Data.Graph.Indexed.Ring.Util
  7 | import Data.SortedMap
  8 | import Derive.Prelude
  9 |
 10 | %language ElabReflection
 11 | %default total
 12 |
 13 | record Profile where
 14 |   constructor P
 15 |   elem      : Elem
 16 |   hcnt      : Nat
 17 |   chrg      : Charge
 18 |   aromBonds : Nat
 19 |   in3ring   : Bool
 20 |   bnds      : Bonds
 21 |
 22 | %runElab derive "Profile" [Show,Eq,Ord]
 23 |
 24 | contribs : SortedMap Profile Double
 25 | contribs =
 26 |   SortedMap.fromList
 27 |     [ (P N 0 0 0 False $ BS 3 0 0, 3.24) -- 1
 28 |     , (P N 0 0 0 False $ BS 1 1 0, 12.36) -- 2
 29 |     , (P N 0 0 0 False $ BS 0 0 1, 23.79) -- 3
 30 |     , (P N 0 0 0 False $ BS 1 2 0, 11.68) -- 4
 31 |     , (P N 0 0 0 False $ BS 0 1 1, 13.6) -- 5
 32 |     , (P N 0 0 0 True $ BS 3 0 0, 3.01) -- 6
 33 |     , (P N 1 0 0 False $ BS 3 0 0, 12.03) -- 7
 34 |     , (P N 1 0 0 True $ BS 3 0 0, 21.94) -- 8
 35 |     , (P N 1 0 0 False $ BS 1 1 0, 23.85) --9
 36 |     , (P N 2 0 0 False $ BS 3 0 0, 26.02) -- 10
 37 |     , (P N 0 1 0 False $ BS 4 0 0, 0.0) --11
 38 |     , (P N 0 1 0 False $ BS 2 1 0, 3.01) --12
 39 |     , (P N 0 1 0 False $ BS 1 0 1, 4.36) --13
 40 |     , (P N 1 1 0 False $ BS 4 0 0, 4.44) --14
 41 |     , (P N 1 1 0 False $ BS 2 1 0, 13.97) --15
 42 |     , (P N 2 1 0 False $ BS 4 0 0, 16.61) --16
 43 |     , (P N 2 1 0 False $ BS 2 1 0, 25.59) --17
 44 |     , (P N 3 1 0 False $ BS 4 0 0, 27.64) --18
 45 |     , (P N 0 0 2 False $ BS 0 0 0, 12.89) --19
 46 |     , (P N 0 0 3 False $ BS 0 0 0, 4.41) --20
 47 |     , (P N 0 0 2 False $ BS 1 0 0, 4.93) --21
 48 |     , (P N 0 0 2 False $ BS 0 1 0, 8.39) --22
 49 |     , (P N 1 0 2 False $ BS 1 0 0, 15.79) --23
 50 |     , (P N 0 1 3 False $ BS 0 0 0, 4.1) --24
 51 |     , (P N 0 1 2 False $ BS 1 0 0, 3.88) --25
 52 |     , (P N 1 1 2 False $ BS 1 0 0, 14.14) --26
 53 |     
 54 |     , (P O 0 0 0 False $ BS 2 0 0, 9.23) --27
 55 |     , (P O 0 0 0 True $ BS 2 0 0, 12.53) --28
 56 |     , (P O 0 0 0 False $ BS 0 1 0, 17.07) --29
 57 |     , (P O 0 (-1) 0 False $ BS  1 0 0, 23.06) --30
 58 |     , (P O 1 0 0 False $ BS 2 0 0, 20.23) --31
 59 |     , (P O 0 0 2 False $ BS 0 0 0, 13.14) --32
 60 |    
 61 |     , (P S 0 0 0 False $ BS 2 0 0, 25.3) --33
 62 |     , (P S 0 0 0 False $ BS 0 1 0, 32.09) --34
 63 |     , (P S 0 0 0 False $ BS 2 1 0, 19.21) --35
 64 |     , (P S 0 0 0 False $ BS 2 2 0, 8.38) --36
 65 |     , (P S 1 0 0 False $ BS 2 0 0, 38.8) --37
 66 |     , (P S 0 0 2 False $ BS 0 0 0, 28.24) --38
 67 |     , (P S 0 0 2 False $ BS 0 1 0, 21.7) --39
 68 |   
 69 |     , (P P 0 0 0 False $ BS 3 0 0, 13.59) --40
 70 |     , (P P 0 0 0 False $ BS 1 1 0, 34.14) --41
 71 |     , (P P 0 0 0 False $ BS 3 1 0, 9.81) --42
 72 |     , (P P 1 0 0 False $ BS 3 1 0, 23.47) --43
 73 |     ]
 74 |
 75 | parameters {0 b,e,p,r,t,ch,l : Type}
 76 |            {auto ce : Cast e Elem}
 77 |            {auto cb : Cast b BondOrder}
 78 |            {k       : Nat}
 79 |            (g       : IGraph k (AromBond b) (Atom e Charge p r HCount t ch l))
 80 |
 81 |   contrib : Elem -> Fin k -> Double
 82 |   contrib el n =
 83 |     let A a ns := adj g n
 84 |         implHs := a.hydrogen
 85 |         hc     := cast implHs.value + countElems g H ns 
 86 |         tm     := inThreeMemberedRing g n
 87 |         pro    := P el hc a.charge (count arom ns) tm (nonAromaticBonds implHs ns)
 88 |      in fromMaybe 0.0 $ lookup pro contribs
 89 |
 90 |   ||| Computes the contribution to the total polar surface area (TPSA)
 91 |   ||| from a single atom.
 92 |   export
 93 |   tpsaContrib : Fin k -> Double
 94 |   tpsaContrib n =
 95 |     case cast @{ce} $ elem $ lab g n of
 96 |       O => contrib O n
 97 |       N => contrib N n
 98 |       S => contrib S n
 99 |       P => contrib P n
100 |       _ => 0.0
101 |
102 |
103 |   ||| Computes the total polar surface area (TPSA) of a molecule.
104 |   export
105 |   tpsa : Double
106 |   tpsa = foldl (\x,v => x + tpsaContrib v) 0.0 (nodes g)
107 |