WebPPL is probably positioned as an educational framework to teach probabilistic programming but I found it has lots of features which makes it ideal for experimentation before moving on to more robust things, like PyMC3 and Pyro. This JavaScript library boasts:

As an example, this is a simple weather model with temperatures depending on a Bernoulli flip whether it’s clouded or sunny:

var weather = function(){
  var maybeCloudy = (flip()? 'cloudy':'sunny');
  var temp = maybeCloudy==='sunny'? gaussian(30,3): gaussian(20,3);  
  var isSunny = maybeCloudy === 'sunny';
  return {sunny: isSunny, temp: temp}
}
 
var d = Infer({}, weather);
viz(d)

It results in a clean line-chart

If you want add a factor (i.e. constraint) you can do so with the ‘condition’ statement:

var weather = function(){
  var maybeCloudy = (flip()? 'cloudy':'sunny');
  var temp = maybeCloudy==='sunny'? gaussian(30,3): gaussian(20,3);
  condition(temp>=24.5 && temp <=25.5)
  var isSunny = maybeCloudy === 'sunny';
  return {sunny: isSunny, temp: temp}
}
 
var d = Infer({}, weather);
viz(d)

corresponding to an observation of a temperature around 25 degrees:

telling you that around this temperature there is an equal chance of having a ‘cloud’ or ‘sunny’ variable. The viz framework auto-detects what works best. If you return a boolean you’ll get a bar-chart which shows even more clearly the result:

var weather = function(){
  var maybeCloudy = (flip()? 'cloudy':'sunny');
  var temp = maybeCloudy==='sunny'? gaussian(30,3): gaussian(20,3);
  condition(temp>=24.5 && temp <=25.5)
  var isSunny = maybeCloudy === 'sunny';
  return   isSunny 
}
 
var d = Infer({}, weather);
viz(d)

Now, let’s do the same with PyMC3. It takes more code and a different way of thinking too:

import pymc3 as pm
import numpy as np

import matplotlib.pyplot as plt
import seaborn as sea
sea.set()
sea.set
%matplotlib inline

weather_model = pm.Model()
with weather_model:
    sunnyTemp = pm.Normal("sunnyTemp", 30,3)
    cloudyTemp = pm.Normal("cloudyTemp", 20,3)
    isCloudy = pm.Bernoulli("isCloudy", p=.5)

    temp1 = pm.math.switch(isCloudy, cloudyTemp, sunnyTemp)
    temp = pm.Normal("temp", temp1, 0.1, observed=[25.])
    step = pm.Slice()
    trace = pm.sample(5000, step=step)
pm.traceplot(trace)

plt.hist(trace['isCloudy'])

In Pyro it amounts to something like

def weather():
    maybeCloudy = 'cloudy' if dist.Bernoulli(.5).sample().item()==1. else 'sunny'
    temp = dist.Normal(30,3) if maybeCloudy == 'sunny' else dist.Normal(20,3)
    isSunny = maybeCloudy=='sunny'
    yn = 1. if isSunny else 0.
    return torch.tensor([yn, temp.sample().item()])
conditioned_scale = pyro.condition(weather, data={"measurement": measurement})

posterior = pyro.infer.Importance(weather, num_samples=100)

marginal = pyro.infer.EmpiricalMarginal(posterior.run())
print(marginal()[1].item())
plt.hist([marginal()[1].item() for _ in range(10000)] ,bins=30)

Although Pyro is the cool new kid in the PPL landscape I found it more abstruse, especially trying to understand how to impose constraints.

Take a look at the wonderful WebPPL tutorials related to multi-agent modeling and how HTML/JS works as a rapid prototyping environment.