Rust for linear algebra and neural networks

From time to time I do a little bit of fooling around with algorithms and AI, I suppose it's a break from what I do in work, which is more about building business applications. I don't have a real goal in mind, I don't think I'm going to stumble about an algorithm for general intelligence or something, but I like to understand how things work at a low level.

So I've writing a bit of Rust code to implement some linear algebra (some matrix and vector operations) and use that to build some small neural networks. You can find the current state of that code here. Of course, if you really want to use production ready code for linear algebra you should check https://www.rustsim.org/ and if neural networks and AI in Rust are what you're after please check https://github.com/AtheMathmo/rusty-machine. Note that rusty-machine is still in a state of flux since it's using its own linear algebra library but there's talk about using nalgebra from rustsim hopefully that will bring together all these initiatives in one unified Rust platform for numeric operations and AI.

I have to say I'm still impressed with Rust. I like that you can do a lot of things with a functional approach, while still being able to update variable values or do other imperative operations when it makes sense. It does feel a lot like "pragmatic Haskell" at times. So far I haven't struggled with the borrow checker, but the code I've written so far is fairly simple and self-contained.

As always, I'm still wondering if there are some projects that it would make sense for me to contribute too. I'm certainly not a AI or math expert, but if you have a Rust project that could do with some contributions, let me know!

Happy Rust Hacking!

Symbolic differentiation to the rescue! ... Or not.

I'm still playing with my LSTM networks, inspired by a few blog posts about their "unreasonable efficiency". I spent a lot of time messing with genetic algorithms, but then I came back to more predictable methods, namely gradient descent. I was trying to optimize the performance of the cost function I use via AD (even asking help on stack overflow) when I stumbled across a couple of blog posts on symbolic differentiation (here and here). The last one, combining automatic and symbolic differentiation, struck a chord. If my differentiation calculation was taking so much time to calculate, could I just not calculate it once with symbolic expressions, then close the resulting expression over my variables (my LSTM network) repeatedly while applying the gradients. I should only pay the price for the derivation once!

So I extended the data type suggested in the blog post to include all operations I was using in my function, manage to sort out all the types and verify via a few tests that it would work. I had great hopes! But as soon as I started testing on a real LSTM, the code just crawled to a halt. I even thought I had some infinite loop, maybe some recursion on the expression, but testing more thoroughly showed that it was the sheer size of the generated expression that was the issue. A LSTM of 2 cells is represented in the cost function used by AD as an array of 44 doubles, so basically for a LSTM of 2 cells, I'll have 44 variables in my gradient expression. My simple test that tries to use a LSTM to generate the string "hello world!" uses 9 cells (9 different characters in the string) , which is 702 variables. Even printing the expression takes forever. Running it through a simplifying step takes also forever. So my idea was not as good as it first looked, but it was fun testing it!

The code for my expressions can be found here, and the code doing the gradient descent via the symbolic differentiation is here. All of that looks probably very naive for you calculus and machine learning experts but hey, I'm a human learning...  If everybody has any idea to speed up my code, I'l happily listen!

Playing with Recurrent Neural Networks in Haskell

Some time ago an interesting article surfaced on Reddit, about using recurrent neural networks to generate plausible looking text. Andrej did a very good job in explaining how they work and some of the techniques and algorithms he used. I thought this was worth some explorations on my own, so I did a bit of research and tried to implement my own networks in Haskell.

I implemented the basics using the hmatrix package to have vectors and matrices. I went to Khan Academy to learn more about these topics because I had stopped math in school before learning matrices. Once I managed to implement hopefully correctly the algorithm I used the simple-genetic-algorithm package (well, my fork, uploaded on hackage) to evolve the network via genetic selection and mutation.

This gave some results, especially when I fixed my mutation code, that in fact was doing nothing, thus emphasizing again that mutation is a critical part of genetic algorithms, and not just crossovers.

Then I went to implement proper learning algorithms, a bit less random than genetic evolution. I was clearly out of my depth there, but learned a lot. To implement gradient descent I actually used the ad library and had to rewrite the algorithms without hmatrix so they would work on simple lists. Given that a couple of weeks again I didn't even understand what AD was, I'm happy I got some stuff to work. I was lucky to find some good code in python, even trying to understand the under-specified RMSProp algorithm.

The end result is not great, though. My network can learn the really complex sentence "Hello world!" and regenerate it after a couple of thousand iterations of the learning algorithm, in around a minute on my machine. I haven't managed to parallelize the treatment yet, and running on all my cores instead of just one make performance a lot worse. On more complex sentences performance becomes really bad, as both the network becomes bigger (more characters to deal with, etc) and the data to evaluate is bigger.

Trying to use a sparse representation for characters using 2 values instead of 1 as in the standard 1-of-k encoding didn't work well. I also realize that tweaks in the cost function had a huge impact on how well and how fast the network can learn. I have started to look at cutting the data in batches but I think I need to make the learning on simple data as fast as possible before moving on.

So, nothing earth shattering, no great progress or insight, but I feel I learned a lot by implementing my own stuff and running into all the issues that are probably well known of people familiar with these subjects. My artificial intelligence is still pretty dumb, but hopefully my own natural intelligence has progressed a bit!

The code is on github of course, so if anybody wants to have a look, I'll gladly take any feedback!!

Genetic evolution of a neural network

In  a previous post I was trying the LambdaNet library for neural networks, training one with a function my own little human brain designed. But can we make the network learn without knowing the actual algorithm?

As a implementation of a genetic algorithm, I use simple-genetic-algorithm, that is, ahem, simple to use compared to some other packages. I've modified it, though, to use MonadRandom instead of having to thread through the RandomGen, so if you want to run the code make sure to take the version from my fork in the MonadRandom branch, until these changes are released by the maintainer on Hackage.

To use a genetic algorithm, we need two things: a representation of the structure we want to evolve that is suitable for the genetic transformations, and the implementation of these transformations (crossover and mutation) themselves.

For the first task, we convert the network to simple vectors:

-- | Convert a network for only vectors structures
toVectors :: (Floating (Vector a), Container Vector a, Floating a) => 
  Network a -> [([Vector a],Vector a)]
toVectors = map (tovs . layerToShowable) . layers
  where
    tovs (ShowableLayer ws bs) = (toRows ws,bs)

-- | Convert back to a network given a network definition and vector structures
fromVectors :: (Floating (Vector a), Container Vector a, Floating a) => 
  LayerDefinition a -> [([Vector a],Vector a)] -> Network a 
fromVectors ld = Network . map ((\sl-> showableToLayer (sl,ld)) . frvs)
  where
    frvs (vs,v)=ShowableLayer (fromRows vs) v

Note that we need a LayerDefinition to rebuild the network from the Vectors. Currently each layer of the network has the same definition and the algorithm does NOT evolve this structure, only the weights and biases.

We're going to keep that information along with some Worlds that we use for fitness testing:

-- | Store the network information as vectors
data NetworkData = NetworkData [([Vector Float],Vector Float)] [World]
  deriving (Show,Read,Eq,Ord)


Then we can implement the Chromosome class. Crossover takes the average of all weights and biases of the two parents, mutation changes one value randomly. Of course other implementations could be found.

-- | Chromosome instance
instance Chromosome NetworkData where
    -- Take average
    crossover (NetworkData xs ws) (NetworkData ys _) =
        return [ NetworkData (Prelude.zipWith zipW xs ys) ws]
        where
          zipW (vs1,v1) (vs2,v2) = (Prelude.zipWith zipW1 vs1 vs2,zipW1 v1 v2)
          zipW1 = V.zipWith (\x y -> (x + y) / 2) 

    -- Mutate one weight randomly
    mutation (NetworkData xs ws) = do
        xs' <- font="" r1="" randomchange="" xs="">
        return $ NetworkData xs' ws
      where
        randomChange _ [] =  return []
        randomChange f xs2 = do
          idx <- -="" 1="" font="" getrandomr="" length="" xs2="">
          mapM (\(i,x)->if i==idx then f x else return x) $ zip [0..] xs2
        r1 (vs,v) = do
          (v3:vs2) <- font="" r2="" randomchange="" v:vs="">
          return (vs2,v3)
        r2 v2 = do
          idx2 <- -="" 1="" font="" getrandomr="" v.length="" v2="">
          dx   <- 20="" font="" getrandomr="" nbsp="">
          return $ v2  V.// [(idx2,dx)]         

    -- calculate fitness
    fitness (NetworkData xs ws) = sum (map (fitWorld xs) ws) / fromIntegral (length ws)

For the fitness function we calculate the fitness for each given world and average it. I'm not trying to be too clever with that code.
For each world we run the food searching algorithm from each corners, and evaluate how far we are from the target, and if we reached it how long it took us. So networks that find the food will always rank higher than the ones who don't, and the quicker among them will rank higher again.

-- | Calculate fitness on a given world   
fitWorld :: [([Vector Float],Vector Float)] -> World -> Double
fitWorld dat w = sum (map fitCorner $ corners $ wSize w) / 4
  where
    fitCorner pos = 
      let network = fromVectors layerDef dat
          poss = algSteps w (neuralAlg w network) 50 pos
          endSmell = currentSmell w $ last poss
          possLength = length poss
      in fitFormula w endSmell possLength (distance pos $ wFood w)
      
-- | Fitness formula
fitFormula :: World -> Int -> Int -> Int -> Double
fitFormula w endSmell possLenth dist = case fromIntegral endSmell / fromIntegral (wSmell w) of
    1 -> 2 + (fromIntegral dist / fromIntegral possLenth)
    n -> n

Then we just need a stop condition: either a maximum number of generations or reaching the maximum possible fitness (shortest path found from all corners)

-- | Maximum for the fitness
maxFit :: Double
maxFit = fitFormula (World undefined undefined 10 undefined) 10 10 10

-- | Stop function
stopf ::  NetworkData -> Int -> IO Bool
stopf nd gen= return $ gen > 300 || fitness nd == maxFit

And code to generate random networks

-- | Build a random network data
buildNetworkData :: (Monad m,RandomGen g) => [World] -> RandT g m NetworkData 
buildNetworkData ws= do
  g <- font="" getsplit="">
  let n = buildNetwork g
  return $ NetworkData (toVectors n) ws

We can then evolve our population of networks, say on two worlds, w and w2:

runGAIO 64 0.1 (buildNetworkData [w,w2]) stopf

 And we get after a couple of minutes a network that can find the food from another point that the tested corners:

Iteration:1

##########

#........#

#........#

#........#

#...X....#

#........#

#........#

#........#

#.......@#

#........#

#........#

##########

Iteration:2

##########

#........#

#........#

#........#

#...X....#

#........#

#........#

#......@.#

#........#

#........#

#........#

##########

Iteration:3

##########

#........#

#........#

#........#

#...X....#

#........#

#.....@..#

#........#

#........#

#........#

#........#

##########

Iteration:4

##########

#........#

#........#

#........#

#...X....#

#....@...#

#........#

#........#

#........#

#........#

#........#

##########

Iteration:5

##########

#........#

#........#

#........#

#...@....#

#........#

#........#

#........#

#........#

#........#

#........#

##########

Yummy!

Interested in all this? I recommend the ai-junkie tutorial, then! My code is of course on github. Feedback and suggestions welcome!

Searching for food using a LambdaNet neural network

From time to time I have a fit and start doing a bit of AI. I saw that a new version of LambdaNet was released so I though I would take it for a spin and try something (a little) bit more complicated that their XOR example.

The problem is simple. In a rectangular world, there is food is one place. The food "smells" and so each position in the world has a smell associated with it, the higher the smell meaning the closer to the food. Can we have a neural network that can navigate to the food?

A few definitions:

-- | Direction to go to
data Direction = TopLeft | Left | BottomLeft | Top | Center | Bottom | TopRight | Right | BottomRight
  deriving (Show,Read,Eq,Ord,Bounded,Enum)

-- | Strength of the smell
type Smell = Int
-- | Input information
type Input = [(Direction,Smell)]
-- | Position in world
type Position = (Int,Int)
-- | Size of the world
type Size = (Int,Int)
-- | Maximum number of steps to take
type Max = Int

-- | Number of directions
dirLength :: Int

dirLength = 1 + fromEnum (maxBound :: Direction)

-- | The world
data World = World
    { wSize   :: Size -- ^ size
    , wSmell  :: Smell -- ^ Smell of the food position
    , wSmells :: DM.Map Position Smell -- ^ All smell strengths by position
    }
    deriving (Show,Read,Eq,Ord)

-- | Function deciding in which direction to move
type StepFunction = Position -> Input -> Direction

Fundamental is the concept of Direction, since we want to move. Basically, when we are in a given position in a world, we can get nine directions and their associated smell (staying in the same place is one position). The function to decide what to do in a given position given all the smells of the neighbouring positions is called StepFunction.

The algorithm is easy to write for a human brain:

-- | The base algorithm: just go toward the highest smell
baseAlg :: StepFunction
baseAlg _ = fst . maximumBy (comparing snd)

Note that we ignore the current position, we only work with the input structure.

On top of that, we need function to build the world with the proper smell indicators, run the algorithm till we find the food, etc. All this code can be found in the GitHub project but is not really critical for our understanding of neural networks. One function of interest is running one step of the algorithm, showing the intermediate structures generated:

-- | Perform one step and return the information generated: direction/smell input, direction output
algStepExplain :: World -> StepFunction -> Position -> (Position,([(Direction,Smell)],Direction))

We get the position back, and the second element of the tuple is the input and the output of the StepFunction.

What we want to do is train a neural network, which should be easy since we have an algorithm we know will work well to find the best position to move to, and then use that network as an implementation of StepFunction.

The hardest in neural network programming is to design the input and output structures, so that they represent adequately the information about your problem in a format that the network can deal with. Here, we have a fixed input size: the smells of the 9 neighbouring positions. The StepFunction returns a Direction, and a Direction is an enum of nine values, so the output of the network could also be 9 values, the highest of these indicating the direction chosen by the network.

The networks in LambdaNet requires Vectors as their input and output data, so lets format the inputs:

-- | Format the inputs suitable for the network
formatInputs ::  World -> [(Direction,Smell)] ->  Vector Float
formatInputs w =   fromList . map (\i-> fromIntegral (snd i) / fromIntegral (wSmell w))    

So an input of 1 means we're on the food itself, and the input value will decrease as we're further from the food, while staying between 0 and 1.

If we have a network, the implementation of StepFunction is straightforward:

-- | Use the network to give the answer 
neuralAlg ::  World -> Network Float -> StepFunction
neuralAlg w n _ is = toEnum $ maxIndex $ predict (formatInputs w is) n 

We format the input, run predict, retrieve the index for the maximum value in the output vector, and use that as the index in the Direction enum. We just need a trained network!

To get that, we generate the training data from a given world. We list all possible positions in the world, calculate the corresponding inputs, run the basic algorithm on the input to get the optimal answer. For the result direction will set the output value to 1, and zero for all the others

-- | Training data: for each position in the world, use the base algorithm to get the training answer
trainData ::  World -> [(Vector Float, Vector Float)]
trainData w = map onePos $ allPositions w
  where
    onePos p = 
      let (_,(is,dir)) = algStepExplain w baseAlg p
          os = map (\(d,_)->if dir==d then 1 else 0) is 
      in (formatInputs w is,fromList os) 

From here, we unimaginatively reuse the LambdaNet tutorial code to build a network...

-- | Create the network
buildNetwork :: RandomGen g => g -> Network Float
buildNetwork g = createNetwork normals g $ replicate 3 $ LayerDefinition sigmoidNeuron dirLength connectFully

And train it:

-- | Train a network on several given worlds
train :: Network Float -> [World] -> Network Float
train n ws = 
  let t = BackpropTrainer (3 :: Float) quadraticCost quadraticCost'
      dat = concatMap trainData ws
  in trainUntilErrorLessThan n t online dat 0.01

What is critical here is that we train the network on several different worlds. I tried training only one world, the resulting network would perform well on worlds of the same size or smaller, but not bigger worlds, because it was too fit for the actual smell values. Training even on only two quite different worlds brought big enhancements in the intelligence of the network, at the code of longer learning time.

Once the network is trained, you can run it on several different worlds and see how it can find the food. There is a simple visualization module that allows you to see clearly the moves, for example:

Iteration 1
##########
#........#
#........#
#........#
#...X....#
#........#
#........#
#........#
#.......@#
#........#
#........#
########## 

(X being the food, @ the current position)

Iteration 3
##########
#........#
#........#
#........#
#...X....#
#........#
#........#
#......@.#
#........#
#........#
#........#
##########

Iteration 6
##########
#........#
#........#
#........#
#...@....#
#........#
#........#
#........#
#........#
#........#
#........#
##########

Yummy!

If you're interested, the full source code with tasty unit tests in on Github.

This is of course very basic, and only begs to be enhanced with more complicated worlds (maybe with walls, several sources of food, places with no smell at all, etc). What do you do when you don't know the best algorithm yourself? Maybe I'll come back for more later to find out!

Circling through ideas

I can't make my mind up about what to do in my free time (I have some, yes yes). I keep circling between the same ideas, thinking "that would be good", then jumping to the next one.

• Keep on working on EclipseFP and other Haskell IDE related ideas. I would like for example to unify the storage of metadata between buildwrapper, scion-browser and eclipsefp, to have one database to would keep library information (definitions, documentations) like scion-browser, AST with types like buildwrapper, usage references like the usage DB in EclipseFP. It could be interesting to try to use that database to drive an IDE and have a clear repository of metadata. But then I'm tired of working on my own on EclipseFP, and when I see that Leksah really has also one active maintainer, I think people are really not interested in advancing Haskell IDEs, they must be happy with Emacs/vi and ghc-mod at a push, so why bother?
• Then I think working on games is fun! I had fun writing Mazes of Monad, and other little games, and I do enjoy playing role playing or adventure games (I can really recommend the last one I've completed, the Longest Journey), so maybe writing a game the reactive way or even for Android like the guys at Keera Studios do would be a good use of my time. Then I remember I'm a programmer that sucks at graphics design and would probably suck as much at game design.
• Games are too trivial, let's do something that will change the world, like work on AI! A few books I've read like Kurzweil's and Hawkins' have been truly inspirational, so maybe I could write some HMM neural temporal gizmo that would become sentient over night!! Then I wake up. People smarter than me and with more time than me are already working on that, so I would not contribute anything anyway. Why don't I help these people by providing a better programming experience, for example a Haskell IDE? Back to idea 1!

Ah well, I'll go back to browsing Reddit!

Fame at last!

I was reading the book "Haskell Data Analysis Cookbook" when suddenly, my name pops up! Funny to see a link to a 7 year old blog entry, who knew I would go down in history for a few line of codes for a perceptron? It's deep in Chapter 7, for those interested. Maybe this is a sign that I should abandon everything else and spend my time on AI, since it's obviously where fame and riches abound! Right...

Hidden Markov Models for Natural Language Tagging, in Haskell

I became intrigued in Hidden Markov Models after reading Kurzweil, who claims they can be used to model the thinking process in the brain. There is much debate about that, but these are interesting AI structures. This page I think has a good introduction.

I was working though the (partial) online book on Natural Language Processing with Haskell, and thought of combining the two. I used Mike Izbicki's hmm library for a one order Hidden Markov Model implementation. Once I initialized the model properly using the training data, I got around 91% accuracy on tagging, which is on par with the rules based approach presented in the nlpwp book.

I used the strategy outline in this paper to deal with unknown words (words in the test set not met in the training set): replace these words with a token that is also used for low frequency words in the training set. So far I've used only one token but I suppose being a bit more fined grained (to distinguish words starting with a capital letter, currency amounts, numbers) will improve results.

Performance is not very good even with some parallelism, so I think I need to spend more time on it, but it's definitely encouraging. It'll be a little bit of time till I have a thinking brain, though, but there is hope!

OpenAlchemist: the start of an AI in Haskell

Lately I've been quite addicted to OpenAlchemist, a Tetris-like game developed by I think a couple of French developers (cocorico!). The code is in C++, there doesn't seem to be a clear API to integrate with the engine, the forums are closed, but I still wanted to write an AI for the game..

So I've written a little something in Haskell and posted it on Github just in case it's of interest to somebody. I basically use the Win32 API to capture the game window and save it as a BMP, then I use the bmp library to open the BMP, and some home grown image recognition to see which shapes are displayed on the screen. Then comes the "easy" part, checking the possible combinations and seeing which one yields the higher score. So far the AI hasn't beaten my own high score, so there's still work to do!

I'm still struggling with shape recognition sometimes, so this is an area that needs an improvement. I probably need to read up a bit on the proper algorithms that people have found to do this kind of things instead of just writing my own, but hey, one needs to have fun!

Evolving a computer with TECS and Genprog

I've started reading The Elements of Computer Systems since it was recommended on Reddit and it's a great read. After the section on boolean gates I decided to have a quick go at build the gates in Haskell. Just implement Nand with pattern matching, and implement all the other gates just in term of functions calls to previously defined gates. The book mentioned "think positive" about implementing Not, so I just T as the second parameter, did't think more about it, even though I wasn't sure having a constant value was OK.



Then I moved on to the other gates. It was interesting to think about boolean logic without Or, for example, after years of only reasoning in terms of And Or and Not. But at some stage I thought "why am I doing this by hand?". I remembered a thread about a new genetic programming library on an Haskell group. Surely generating the gates would be simple enough?
So I got genprog from hackage, and adapted the example for my needs. I started with a very simple system to determine the implementation of Not:



And of course, the result came quickly enough: not a=nand(a,a). No need for a Const at all! Evolution has beaten me! (OK, I'm only a not-so-intelligent designer, I suppose...)
Now, I'm going to see if this approach can take me all the way to the arithmetic unit...

Haskell Neural Network: plugging a space leak

First, the good news: last week when I posted about my little digit recognition program I had made a mistake in my test code: when I give a rounded eight to my network, it does recognize it, even though it's been trained only with a square eight! Yoohoo!


Also, a few people pointed to online database of handwritten digits so that I could get more data to train my network. I dully downloaded and started working with some files. Soon I knew I had a problem. Not only was the learning slow, but the memory usage of the program showed exponential growth, like using 400MB after a few iterations. Ha ha, I thought, so it's me against lazy evaluation again!


Well, I was right. A couple of strictness annotation in the neuron output calculation function helped performance, nearly making the code twice as fast, but didn't help the memory consumption. Going mad with ! (BangPatterns is your friend) didn't help any more. WTF?


Then I started to read a bit more about seq and Weak Head Normal Form. So when I put a BangPattern in front of a tuple, it evaluates the tuple but not what's inside it. And, lo and behold, if at each learning iteration I traced the network to the standard output, memory usage stayed constant! So I needed to tell Haskell to fully evaluate the network each time, so as not to keep thunks full of arithmetic operations hanging round. A distant memory of a package called deepseq provided me the final help I needed. I didn't use the package itself because I wanted to understand really what's going on, so I just roll out a few lines of my own:



And by calling deepSeq on my network in each iteration, my program runs in constant space! And I looked at the flat line of memory usage in my task manager, and I was pleased. Now, I can try to improve the speed of the dang thing, but that's another set of problems...

Digit recognition with a neural network. First attempt!

I bought myself "On Intelligence", and I thought as the first step on the long road to building an intelligent machine I would go back to neural networks. I have my own little library that seems to work (kind of), so I decided I was going to try to write a little number recognition program. I was inspired of course by ai-junkie, but it was a bit short on details...
So here is the problem: I have a 5x5 grid of LEDs that can be on or off, and I want to recognize the digit that they show. 
The input will be a list of 25 numbers, either 0 (off) or 1 (on)
The output will be a list of 10 numbers, either 0 or 1. The first position of a 1 indicates the result (so if the first element is a 1, it's number 0, etc).
I arbitrarily decide to set the number of hidden neurons to 50, twice the input neurons. So my topology is 25-50-10. 
I train the network with only one version of each digit. For example 3 would be:

*****
    *
 ****
    *
*****

Wrote a few lines of Haskell code to read the training sets (moving from the representation with *, spaces and new lines to the list of 1 and 0, so that 3 becomes: [1.0,1.0,1.0,1.0,1.0,0.0,0.0,0.0,0.0,1.0,0.0,1.0,1.0,1.0,1.0,0.0,0.0,0.0,0.0,1.0,1.0,1.0,1.0,1.0,1.0]).
Then I train the network. It learns after 200 iterations, and successfully recognizes the training sets (ok, so things are working as they should be).
Then the real test: can my network recognize small variants of the training set?
Lets's try with:

****
    *
 ***
    *
****

(a slightly rounded 3): the network answers 3!!
Works also for a slightly rounded 6. But then:

 ***
*   *
 ***
*   *
 ***

(rounded 8) gives me 6. Oh dear... I could probably add more cases in the training sets and hopefully resolve the ambiguities, but there are bigger issues (I know I know, I'm very candid, but it's one thing being told something in a book, and another to build something to actually show it to you). For example:

***
* *
***

Is a small zero that doesn't take the full 5 LED width. My network recognizes 0 when it fills the full square, not when it's smaller (tells me it's a 5).
And of course, 1 is a disaster. The training "1" is a vertical row of LEDs on the left. A vertical row of LEDs on the right gives me 9, usually (and I suppose, 9 has all the LEDs on the right on).
So not only is basic training not sufficient to detect simple variations, but also the fact that we deal with absolute positioning of leds mean than scaling and simple side translations are not supported. So I suppose the next step is to not work with absolute positions, but relative positions of "on" LEDs. Another day's work...

Found my neural network bug (three years later)

After the university of Waterloo AI Challenge I went back to a bit of AI programming in Haskell. Three years ago I had written some code to implement a neural network, and it didn't seem to learn as fast as the book said it should. I had not worked on that any more, but I decided to look at it again. And I found a bug!! I mixed up the deltas to pass on from one run to the next. Changing two lines of code makes my network learning much faster.

So now onto something a bit more meaty than exclusive or: very simple character recognition. I want to be able to recognize numbers written as a digital alarm clock would display them (at least mine): you have 7 lights, four vertical and three horizontal, and various patterns make the numbers. If all lights are on, it's 8, if only the middle horizontal light is off we have zero, etc...

I just test things with a network that has 7 input neurons (one for each light), 10 output neurons (one for each number) and 8 hidden neurons (the average between input and output)

g<-getStdGen

n<-network 7 8 10

let trainingSet=[

([1,1,1,-1,1,1,1],[1,0,0,0,0,0,0,0,0,0]), -- 0

([-1,-1,1,-1,-1,1,-1],[0,1,0,0,0,0,0,0,0,0]), -- 1

([1,-1,1,1,1,-1,1],[0,0,1,0,0,0,0,0,0,0]), -- 2

([1,-1,1,1,-1,1,1],[0,0,0,1,0,0,0,0,0,0]), -- 3

([-1,1,1,1,-1,1,-1],[0,0,0,0,1,0,0,0,0,0]), -- 4

([1,1,-1,1,-1,1,1],[0,0,0,0,0,1,0,0,0,0]), -- 5

([1,1,-1,1,1,1,1],[0,0,0,0,0,0,1,0,0,0]), -- 6

([1,1,1,-1,-1,1,-1],[0,0,0,0,0,0,0,1,0,0]), -- 7

([1,1,1,1,1,1,1],[0,0,0,0,0,0,0,0,1,0]), -- 8

([1,1,1,1,-1,1,1],[0,0,0,0,0,0,0,0,0,1]) -- 9

]

let (LNS(n',err,ds),e) = run (initialState n) defaultLearningRate trainingSet (1,1000)

And after something like 150 iterations I get a network that has a low error rate. A little helper method to massage the results:

runNumber n inputs=do

let

(_, output)=exec n inputs

normalized=map (\v->if (1-v)>0.5 then 0 else 1) output

return $ elemIndex 1 normalized

Et voila!

Main> runNumber nn [1,-1,1,1,-1,1,1]

Just 3

Now, things are starting to get interesting! Hopefully I won't wait another three years before trying maybe more complex character recognition.

Waterloo Contest and Human Intelligence

I've entered the Waterloo university Google AI Challenge, hiding under a cunning nickname, since I saw that you could develop in Haskell, and a few others Haskellers are on it too! At the moment, I'm cruising at around the 110th position, which is not too bad given that I only worked a few hours at night on it (I have a day job, you know...) and that I implemented everything myself since the starter pack wasn't ready when I started working. I've played a bit with the standard strategies, implemented some flood fill, reuse some A* implementation I have.
All of this is interesting, of course, but there is no AI as such. The only intelligence in my code is the one I put in, and the strategies I devised. So I'm thinking of playing around with genetic programming and evolution techniques to try to get my code to find the right strategies at the right time. I have the feeling the trick will be to design something that is not too restricted and is truly able to evolve into some more advanced and complex than what I've designed myself. But hey, I'll learn, even if my program doesn't (-:

Mazes of Monad: 1.0.5 and (maybe) final version uploaded

I've uploaded yesterday on Hackage the 1.0.5 version of my Role Playing Game called Mazes of Monad. I think this will be the last upload unless some big bugs are found. But I think I'll stop here in terms of functionality. It was fun, I learned a lot about Haskell, both on programming with the language, and on how to design Haskell code. Working on a full project over a long period of time forced me to actually understand code I had written previously and make it evolve instead of just having a clean slate.

I also realized than while developing a game is a lot of fun, the actual programming is the, ahem, easy part. I find that the actual game design (how to get usable and useful features) and the tweaking of the difficulty, for example, are the most challenging. For a while my test character was quite a strong warrior and I thought the game was easy, so I changed some settings to make it harder, then I realized that I had made it quite hard to start with a poor naked wizard that knew only one spell... I can see that balancing all the different game settings and getting the right formulas is critical, and that only comes from testing and collecting feedback...

Anyway, I have ideas for other games, probably the next one will be web based to have some kind of GUI and be multi player, maybe using happstack as the server (or should I just use a little home grown server in Haskell, for fun?). I'll keep on posting about that!

Haskell RPG Game uploaded on Hackage!

I've worked on it on and off for months, I've rewritten parts of it entirely, it gave me the opportunity to knock together a real, complete application in Haskell, and it's finally done! I've uploaded my console based RPG game (read: nethack clone) to Hackage (get it here). You create a character and walk through mazes, trading, stealing, casting spells and fighting monsters.
You can look at the source code (no documentation, sorry) and see what I've used. I've used a monad to abstract away the state and the random generator (so I can plug a non random generator for testing, one that for example will guarantee that the next throw of the D20 die will give 8). I have developped my own little framework for console handling (executing actions from what the user typed in and displaying results) maybe I should look at the existing libraries and try to plug one in instead.
What I found was the most awkward was the record syntax for nested structures, I'll have to investigate to find an easier way to read and update my structures.
There are a few HUnit tests to verify that the most basic functionality works, so hopefully there is no huge bug.
It was fun to play around with Haskell and with game concepts. I hope I'll have time to work on even cooler stuff like adding AI for the monsters (so they can decide when to run away, etc...).
Any feedback on the game mechanism or on the code itself is welcome!

A very dumb neural network in Haskell

Continuing my implementations of Michael Negnevitsky book on Artifical Intelligence in Haskell, I have tried to code a three layer neural network to implement exclusive or (Intelligence indeed...).

A neuron is made of its input weight and its threshold, and its output is calculated thus (apologies to all of you who gave me some tips about parameter order and newtype declaration, I haven't changed the code yet to follow your good advice!):

type Neuron= ([Value],Value)
type Value = Float

neuronOutput :: [Value] -> Neuron -> Value
neuronOutput inputs (weights,threshold) =
    let
        tot=(foldl (+) 0 (zipWith (*) inputs weights)) - threshold
    in
        1 / (1 + ((2.7182) ** (-tot)))

A three layer network whose input layer does nothing but fire its input to the midle layer can be modelled by to lists of neurons, and I use a TestCase type to keep track of inputs and expected outputs

       
type Network = ([Neuron],[Neuron])
type TestCase = ([Value],[Value])

Running the whole network with some input values is easy:

exec ::  Network -> [Value] -> ([Value],[Value])
exec (hidden,output) inputs =
    let
        hiddenOutputs = map (neuronOutput inputs) hidden
        outputOutputs = map (neuronOutput hiddenOutputs) output
    in
        (hiddenOutputs,outputOutputs)

The meat is all in one function, that calculates the output of the network for a given test and adapts the weights accordingly, first for the output layer then for the hidden layer:

defaultLearningRate::Value
defaultLearningRate=0.1

step ::  (Network,Value) -> TestCase -> (Network,Value)
step ((hidden,output),err) (inputs,expected)=
    let
        (hiddenOutputs,outputOutputs)=exec (hidden,output) inputs
        errorOutputs= zipWith (-) expected outputOutputs
        gradientOutputs=zipWith (\o d->o*(1-o)*d) outputOutputs errorOutputs
        newOutputs=zipWith (\(ws,t) g -> ( zipWith (\h w->w+ (defaultLearningRate*h*g)) hiddenOutputs ws ,t+(defaultLearningRate*g*(-1)))) output gradientOutputs
        weightsByHidden= transpose (map fst output)
        gradientHidden =zipWith (\o ws-> o * (1-o) * (foldl (+) 0 (zipWith (*) ws gradientOutputs))) hiddenOutputs weightsByHidden
        newHidden=zipWith (\(ws,t) g -> ( zipWith (\h w->w+ (defaultLearningRate*h*g)) inputs ws ,t+(defaultLearningRate*g*(-1)))) hidden gradientHidden
        newErr=err + (foldl (+) 0 (map (^2) errorOutputs))
    in ((newHidden,newOutputs),newErr)
 

The result is a new adapted network and the sum of squared errors.

An epoch is a set of tests:

epoch :: Network -> [TestCase] -> (Network,Value)
epoch network= foldl step (network,0)

We can then repeat epoch till the sum of squared errors for an epoch is less than the threshold for convergence:

errorConvergence::Value
errorConvergence=0.001

run ::     Network -> [TestCase] -> Int -> (Network,Int,Value)
run network allInputs epochNb=
    let (newNetwork,delta) = epoch network allInputs
    in (?) (delta <= errorConvergence || epochNb>225)
            (newNetwork, epochNb,delta)
            (run newNetwork allInputs (epochNb+1))
   
(?) :: Bool -> a -> a -> a
(?) True a _  = a
(?) False _ a = a

(Is there a standard function for (?)? I got tired of if then else...)

Then I test my code with the same start network as in the book:

test = do
    let n = ([([0.5,0.4],0.8),([0.9,1.0],-0.1)],[([-1.2,1.1],0.3)])
    let (n',e,err) = run n [([1,1],[0]),([0,1],[1]),([1,0],[1]),([0,0],[0])] 1
    print (n',e,err)
    return ()

And I get a mixed feeling about the result: the network does converge, and then it does work, but after something like 54000 epochs, while in the book it supposedly take 224! The book gives the result for the first test and I'm getting exactly the same thing, but it looks like my network learns a lot slower than it should. I would suppose I have a subtle error somewhere, but for the life of me I can't find it... If it jumps at anybody, please let me know. My hope now is that when refactoring the code and following the books advice on improving the performance I will find the bug!!

I have tried with random networks:

network :: Int -> Int -> Int-> IO Network
network inputNb hiddenNb outputNb= do
    a<-randomNeurons inputNb hiddenNb
    b<-randomNeurons hiddenNb outputNb
    return (a,b)
   
randomNeurons:: Int -> Int -> IO [Neuron]
randomNeurons nbInput nbNeurons= replicateM nbNeurons (randomNeuron nbInput)
   
randomNeuron:: Int -> IO Neuron
randomNeuron nb= do
    w<-replicateM nb (randomWeight nb)
    t<-randomWeight nb
    return (w,t)
   
randomWeight:: Int -> IO Value
randomWeight nbInput= do
    let interval= randomR ((-2.4)/(fromIntegral nbInput),2.4/(fromIntegral nbInput))
    getStdRandom interval

But I don't get better results...

A perceptron in Haskell

After reading Michael Negnevitsky book on Artifical Intelligence, I started playing with some of the algorithms he gives in Haskell. So my first endeavour is writing a simple perceptron (a one neuron neural network) to perform AND or OR logical operations. Next I'll start a full multi layer network to do much more complex operations like exclusive OR (-:.
Once again I was struck with the simplicity of the Haskell code necessary to implement such functionality (and I'm no Haskell expert yet). You think about how you going to implement it, the type checker tells you when you're doing something stupid, and by the time everything is written and compiles it kinda works.
The core method is calculating the output of the neuron, given the inputs, the weight of each input and the threshold value for that neuron:

neuronOutput :: [Float] -> [Float] -> Float -> Float
neuronOutput inputs weights threshold =
    let
        tot=(foldl (+) 0 (zipWith (*) inputs weights)) - threshold
    in
        if tot >= 0
            then 1
            else 0

Using zipWith we multiply each input with its weigth, calculate the sum with fold and remove the threshold. A negative result yields 0, else it yields 1.

Adjusting the weights to adapt to the error of a test case is no more difficult:

defaultLearningRate::Float
defaultLearningRate=0.1

adjustWeights :: [Float] -> [Float] -> Float -> Float -> [Float]
adjustWeights inputs origWeights expected actual =
    let
        e=expected - actual
        delta (i,w) = w + (defaultLearningRate * i * e)
    in
        map delta (zip inputs origWeights)

We modify the weight by applying to each weight a delta that is dependant on the input weight, the error ratio and the learning rat

A single step is made of calculating the output of the neuron and adapting the weights, given the expected value:

defaultThreshold::Float
defaultThreshold=0.2

step :: [Float] -> [Float] -> Float -> [Float]
step inputs weights expected=
    let
        o = neuronOutput inputs weights defaultThreshold
    in
        adjustWeights inputs weights expected o   

Then an epoch applies the step function to each different input sets. The hardest part is calculating the average of all deltas to find if the weights have changed at all during the evaluation of all the test cases

epoch :: [([Float],Float)] -> [Float] -> ([Float],Float)
epoch allInputs weights=
    let
        f w (inputs,expected) = step inputs w expected
        newWeights=foldl f weights allInputs
        delta= (foldl (+) 0 (map abs (zipWith (-) newWeights weights))) / (fromIntegral $ length weights)
    in (newWeights,delta)

A training run is only launching the epoch function repeatedly till the delta is zero:

run ::     [([Float],Float)] -> [Float] -> Int -> ([Float],Int)
run allInputs weights epochNb=
    let (newWeights,delta) = epoch allInputs weights
    in if delta == 0
        then
            (newWeights, epochNb)
        else
            run allInputs newWeights (epochNb+1)

We generate the initial weights in the IO monad, passing the number of weights needed:

initialWeights :: Int -> IO [Float]
initialWeights nb = do
    let interval= randomR (-0.5,0.5)
    (replicateM nb (getStdRandom interval))

Putting it all together:

test :: IO()
test = do
    w<-initialWeights 2
    let (ws,i)=run [([0,0],0),([0,1],0),([1,0],0),([1,1],1)] w 1
    print (ws,i)
    return ()

test find out the proper weights to implement AND and print them out, with the number of epochs needed to reach the result. You can then call neuronOutput with the resulting weights and your input: you have a AND logical gate!

As usual, any pointer to how to improve the code is welcome!