added automimicry-effect

app/Main.hs

@@ -25,11 +25,6 @@ pps = Enzyme "PPS" [(Substrate Phosphor,1)] ((Substrate Phosphor,-1),(Produced P
 fpps :: Enzyme -- PP -> FPP
 fpps = makeSimpleEnzyme (Produced PP) (Produced FPP)
 
--- Predator
-
-greenfly :: Predator         -- 20% of plants die to greenfly, but the fly is
-greenfly = Predator [] 0.2   -- killed by any toxic Component
-
 -- Environment
 
 exampleEnvironment :: Int -> [Enzyme] -> [(Predator,Probability)] -> [(Compound,Amount)] -> Environment
@@ -44,6 +39,10 @@ exampleEnvironment addedC es pred tox =
     , maxCompound = maxCompoundWithoutGeneric + addedC
     , toxicCompounds = tox --[(Produced FPP,0.1)] ++ tox
     , possibleEnzymes = es -- [pps,fpps] ++ es
+    , settings = Settings { automimicry = True
+                          , predatorsRandom = False
+                          , numPlants = 150
+                          }
     }
 
 -- Plants
@@ -76,25 +75,37 @@ examplePlants = (\g -> Plant g defaultAbsorption) <$> genomes
 -- ----------------------
 
 loop :: Int -> [Plant] -> Environment -> IO ()
-loop loopAmount = loop' loopAmount 0
+loop loopAmount ps env = loop' loopAmount 0 ps env
 
   where
+    -- cache enzyme colorful-strings
+    stringe :: [(Enzyme, String)]
+    stringe = (\e -> case Data.List.find (\(t,_) -> (t==) . fst . snd . synthesis $ e) toxins of
+               Just (_,toxicity) -> (e,"\ESC[38;5;" ++ show (16 + 36*5 + 6*floor (5*(1-toxicity)) + 0) ++ "m" -- yellow -> red rainbow for tocixity 0 -> 1
+                                       ++ padded 50 (show (enzymeName e)) ++ "\ESC[0m")
+               Nothing           -> (e, padded 50 (show (enzymeName e)))
+             ) <$> possibleEnzymes env
+    toxins :: [(Compound, Amount)]
+    toxins = toxicCompounds env
+    padded i str = take i $ str ++ repeat ' '
+    printEvery = 1
     loop' :: Int -> Int -> [Plant] -> Environment -> IO ()
     loop' loopAmount curLoop plants e = unless (loopAmount == curLoop) $ do
+      when (curLoop `mod` printEvery == 0) $ do
         putStr "\ESC[2J\ESC[H"
         printEnvironment e
         putStrLn ""
         putStrLn $ "Generation " ++ show curLoop ++ " of " ++ show loopAmount ++ ":"
       newPlants <- flip runReaderT e $ do
-          ! fs <- fitness plants
+          ! fs <- fmap (+0.01) <$> fitness plants -- fitness should be at least 0.01 for mating to work
           let fps = zip plants fs -- gives us plants & their fitness in a tuple
               sumFitness = sum fs
-          pe <- asks possibleEnzymes
+          when (curLoop `mod` printEvery == 0) $ do
-          tc <- asks toxicCompounds
+            liftIO $ printPopulation stringe fps
-          liftIO $ printPopulation tc pe fps
             liftIO $ hFlush stdout
-          -- generate 100 new plants.
+          -- generate x new plants.
-          sequence . flip fmap [1..100] $ \_ -> do
+          np <- asks (numPlants . settings)
+          sequence . flip fmap [1..np] $ \_ -> do
                 parent' <- liftIO $ randomRIO (0,sumFitness)
                 let
                     -- if we only have one parent in our list, take it.
@@ -107,21 +118,22 @@ loop loopAmount = loop' loopAmount 0
                     parent = findParent parent' fps
                 haploMate parent
       hFlush stdout
-      threadDelay $ 100*1000 -- sleep 100ms
+      when (curLoop `mod` printEvery == 0) $ do
+        threadDelay $ 100*1000 -- sleep x*1000ns (=x ~ ms)
       loop' loopAmount (curLoop+1) newPlants e
 
 main :: IO ()
 main = do
   hSetBuffering stdin NoBuffering
-  hSetBuffering stdout NoBuffering
+  --hSetBuffering stdout NoBuffering
-  randomCompounds <- makeHead (Substrate Photosynthesis) <$> generateTreeFromList 50 (toEnum <$> [(maxCompoundWithoutGeneric+1)..] :: [Compound]) -- generate roughly x compounds
+  randomCompounds <- makeHead (Substrate Photosynthesis) <$> generateTreeFromList 30 (toEnum <$> [(maxCompoundWithoutGeneric+1)..] :: [Compound]) -- generate roughly x compounds
   ds <- randoms <$> newStdGen
   probs <- randomRs (0.2,0.7) <$> newStdGen
-  let emptyPlants      = replicate 50 emptyPlant
+  let poisonedTree     = poisonTree ds randomCompounds
-      poisonedTree     = poisonTree ds randomCompounds
+      poisonCompounds  = foldMap (\(a,b) -> [(b,a) | a > 0.2]) poisonedTree
-      poisonCompounds  = foldMap (\(a,b) -> [(b,a) | a > 0.5]) poisonedTree
+  predators <- generatePredators 0.5 poisonedTree
-  predators <- generatePredators 0.8 poisonedTree
   let env              = exampleEnvironment (getTreeSize randomCompounds) (generateEnzymeFromTree randomCompounds) (zip predators probs) poisonCompounds
+      emptyPlants      = replicate (numPlants . settings $ env) emptyPlant
   printEnvironment env
   writeFile "poison.twopi" $ generateDotFromPoisonTree "poison" 0.5 poisonedTree
   putStr "\ESC[?1049h"
@@ -130,14 +142,6 @@ main = do
   _ <- getChar
   putStr "\ESC[?1049l"
 
-  -- fitness <- runReaderT (sequence $ (\a -> do p <- absorbNutrients a >>= produceCompounds a; (,,) a p <$> deterPredators p) <$> emptyPlants) exampleEnvironment
-  -- mapM_ (printf "%15.15s, " . show . toEnum @Compound) [0..maxCompoundWithoutGeneric]
-  -- putStrLn "Fitness"
-  -- forM_ fitness $ \(p, c, f) -> do
-  --   mapM_ (printf "%15.2f, ") (toList c)
-  --   printf "%15.2f" f
-  --   putStr "\n"
-
 generatePredators :: Double -> EnzymeTree s (Double,Compound) -> IO [Predator]
 generatePredators threshold t = do
     ps <- mapM generatePredators' $ getSubTrees t
@@ -148,13 +152,13 @@ generatePredators threshold t = do
           let comps = foldMap (\(a,b) -> [(a,b) | a > threshold]) t
           amount <- randomRIO (0,length comps + 1) :: IO Int
           forM [1..amount] $ \_ -> do
-              impact <- randomRIO (0.2,0.7)
+              impact <- randomRIO (0.1,0.2)
               rands <- randoms <$> newStdGen
               let unresists = foldMap (\((a,b),r) -> [b | r*2 < a]) $ zip comps rands
-              return $ Predator unresists impact
+              return $ Predator unresists impact 1
 
 printEnvironment :: Environment -> IO ()
-printEnvironment (Environment soil pred metaIter maxComp toxic possEnz) =
+printEnvironment (Environment soil pred metaIter maxComp toxic possEnz settings) =
   do
     putStrLn "Environment:"
     putStrLn $ "Soil:      " ++ show soil
@@ -162,18 +166,16 @@ printEnvironment (Environment soil pred metaIter maxComp toxic possEnz) =
     putStrLn $ "PSM Iters: " ++ show metaIter
     putStrLn $ "Compounds: " ++ show ((toEnum <$> [0..maxComp]) :: [Compound])
     putStrLn $ "Toxic:     " ++ show toxic
+    putStrLn $ "Settings:  " ++ show settings
 
-printPopulation :: [(Compound, Amount)] -> [Enzyme] -> [(Plant,Double)] -> IO ()
+printPopulation :: [(Enzyme,String)] -> [(Plant,Double)] -> IO ()
-printPopulation toxins es ps = do
+printPopulation es ps = do
   let padded i str = take i $ str ++ repeat ' '
   putStr $ padded 50 "Population:"
   forM_ ps $ \(_,f) -> putStr (printColor f '█')
   putStrLn colorOff
-  forM_ es $ \e -> do
+  forM_ es $ \(e,s) -> do
-    putStr $ case Data.List.find (\(t,_) -> (t==) . fst . snd . synthesis $ e) toxins of
+    putStr s
-               Just (_,toxicity) -> "\ESC[38;5;" ++ show (16 + 36*5 + 6*floor (5*(1-toxicity)) + 0) ++ "m" -- yellow -> red rainbow for tocixity 0 -> 1
-                                 ++ padded 50 (show (enzymeName e)) ++ "\ESC[0m"
-               Nothing           -> padded 50 (show (enzymeName e))
     forM_ ps $ \(Plant g _,_) -> do
       let curE = sum $ map (\(_,q,a) -> fromIntegral q*a)
                      . filter (\(e',_,_) -> e == e')

src/Environment.hs

@@ -84,9 +84,19 @@ data Predator = Predator { irresistance    :: [Compound]
                          , fitnessImpact :: Amount
                            -- ^ impact on the fitness of a plant
                            --   (~ agressiveness of the herbivore)
+                         , numAttacks :: Amount
+                           -- ^ Avarage number of attacks in a generation of appearance
+                           --   (~ mean of poisson-distribution)
                          } deriving (Show, Eq)
 
--- The environment itself is just the soil and the predators. Extensions would be possible.
+-- | Settings to enable/disable parts of the simulation
+data Settings = Settings { automimicry :: Bool     -- ^ do we have automimicry-protection?
+                         , predatorsRandom :: Bool -- ^ do predators always appear or according to their random distribution?
+                         , numPlants :: Int        -- ^ number of plants in starting population
+                         }
+                         deriving (Show, Eq)
+
+-- | The environment itself.
 
 data Environment =
      Environment
@@ -108,6 +118,7 @@ data Environment =
      --   Kills 100% of Predators above Amount.
    , possibleEnzymes :: [Enzyme]
      -- ^ All enzymes that can be created by genetic manipulation in this setting.
+   , settings :: Settings
    } deriving (Show, Eq)
 
 -- helper function. Allows for [0..maxCompoundWithoutGeneric] :: [Compound] with all non-generic Compounds
@@ -151,23 +162,28 @@ fitness ps = do
   products <- sequenceA $ zipWith produceCompounds ps nutrients -- produce compounds
   ds <- liftIO $ randoms <$> newStdGen
   preds <- asks predators
+  randPred <- asks (predatorsRandom . settings)
   let
-     appearingPredators = fmap fst . filter (\((_,p),r) -> p > r) $ zip preds ds -- assign one probability to each predator, filter those who appear, throw random data away again.
+     appearingPredators = if randPred then
-                                                                              -- appearingPredators is now a sublist of ps.
+                            fmap (fst . fst) . filter (\((_,p),r) -> p > r) $ zip preds ds -- assign one probability to each predator, filter those who appear, throw random data away again.
-  survivalRate <- mapM (deterPredators products preds) products  -- defeat predators with produced compounds
+                                                                                           -- appearingPredators is now a sublist of preds without the probability.
+                          else
+                            fst <$> preds -- else just forget about probabilities
+  automimicry <- asks (automimicry . settings)
+  popDefense <- if automimicry then
+                   forM appearingPredators $ \p -> do
+                      as <- mapM (deterPredator p) products      -- how good can an individual deter p
+                      return $ sum as / fromIntegral (length as) -- how good can the population deter p on average
+                else
+                  return $ repeat 1
+  survivalRate <- mapM (deterPredators (zip appearingPredators popDefense)) products  -- defeat predators with produced compounds
   let sumEnzymes = sum . fmap (\(_,q,a) -> fromIntegral q*a) . genome <$>  ps -- amount of enzymes * activation = resources "wasted"
       staticCostOfEnzymes = (\x -> 1 - 0.01*x) <$> sumEnzymes -- static cost of creating enzymes
-      -- primaryEnzymes = filter (\(e,_,_) -> case (fst.fst.synthesis) e of -- select enzymes which use substrate
-      --                                        Substrate _ -> True
-      --                                        otherwise   -> False)
-      --                         (genome p)
   nutrientsAvailable <- fmap snd <$> asks soil
   let nutrientsLeft = (\p -> [p ! i | i <- [0..fromEnum (maxBound :: Nutrient)]]) <$> products
       nutrientRatio = minimum . zipWith (flip (/)) nutrientsAvailable <$> nutrientsLeft
       costOfEnzymes = max 0 <$> zipWith (\s n -> s-n*0.1) staticCostOfEnzymes nutrientRatio -- cost to keep enzymes are static costs + amount of nutrient sucked out of the primary cycle
   return $ zipWith (*) survivalRate costOfEnzymes
-  -- can also be written as, but above is more clear.
-  -- fitness p = absorbNutrients p >>= produceCompounds p >>= deterPredators
 
 produceCompounds :: Plant -> [(Nutrient, Amount)] -> World (Vector Amount)
 produceCompounds (Plant genes _) substrate = do
@@ -186,19 +202,23 @@ produceCompounds (Plant genes _) substrate = do
     -- faster, because no inversions and optimized eig.
   return final
 
+-- Automimicry: see https://www.ncbi.nlm.nih.gov/pmc/articles/PMC2275178/#__sec2title Formula 2.1
+-- Note: F(D) is "costOfEnzymes", but in 'fitness' we multiply "costOfEnzymes" already,
+-- so F(D) is omitted
+-- A(d_hat) is ahat * numAttacks p, because ahat is only deterrence of the population
+-- and does not incorporate the number of attacks, which A(d_hat) in the paper does.
+deterPredators :: [(Predator, Double)] -> Vector Amount -> World Probability
+deterPredators appearingPredators compounds = do
+  deters <- forM appearingPredators $ \(p,ahat) -> do
+    myDeter <- deterPredator p compounds
+    return $ exp $ negate $ numAttacks p * ahat * myDeter -- exp due to assumption that number of attacks are poisson-distributed.
+  return $ product deters
 
--- TODO:
+
--- - dampen full-force due to auto-mimicry-effect. => Fitness would not depend on single plant.
+deterPredator :: Predator -> Vector Amount -> World Double
-deterPredators :: [Vector Amount] -> [(Predator,Amount)] -> Vector Amount -> World Probability
+deterPredator p comps = do
-deterPredators others appearingPredators cs = do
+  toxins <- asks toxicCompounds
-  -- ps <- asks predators
+  return $ product [1 - min 1 (comps ! fromEnum t * l) | (t,l) <- toxins, t `elem` irresistance p]
-  ts <- asks toxicCompounds
-  let
-     deter :: Predator -> Double
-         -- multiply   (toxicity of t with 100% effectiveness at l| for all toxins t; and t in p's irresistance-list)
-     deter p = product [1 - min 1 (cs ! fromEnum t / l)           | (t,l) <- ts,       t `elem` irresistance p]
-         -- multiply  (probability of occurence * intensity of destruction / probability to deter predator | for all predators)
-  return $ product [min 1 ((1-prob) * fitnessImpact p / deter p) | (p,prob) <- appearingPredators]
 
 -- Mating & Creation of diversity
 -- ------------------------------