title: Fitting a Python Scipy truncnorm model to observed values and then sampling
body: <p>I am struggling with the Scipy truncnorm fit method and I would appreciate help so that the fitted parameter coefficients are consistent with the observed data.</p>
  
  <p>As an example, I have created a small sample from the right hand tail of a N(0,1) distribution (where the observations are larger than 2 stdev) and have thrown in a few outliers.</p>
  
  <pre><code>import numpy as np
  import pandas as pd
  import matplotlib.pyplot as plt
  from scipy.stats import truncnorm
  
  values = np.array([2.01,  2.06,  2.71,  2.31,  2.58,  2.17,  3.03,  2.24,  2.12,
                     2.72,  2.46,  2.66,  2.49,  3.41,  2.46,  2.12,  2.12,  2.65,
                     2.32,  2.49,  5.15,  2.62,  2.48,  2.27,  2.05])
  pd.Series(values).describe()
  </code></pre>
  
  <p>This then produces the following summary statistics:</p>
  
  <pre><code>count    25.000
  mean      2.548
  std       0.633
  min       2.01
  25%       2.17
  50%       2.46
  75%       2.65
  max       5.15
  </code></pre>
  
  <p>To illustrate the problem I am having with the scipy fit method and to better understand the truncnorm implementation, I have built the following intuitive models from inspecting the above summary statistics and sampled histograms to the observed values (see charts below). What I am struggling with is why the fit method gives such poor results when I attempt to sample using the estimated parameters? In the event I am not using the fit results correctly or making some other mistake, I would appreciate help with the transformations?</p>
  
  <p>The code to build up these examples:</p>
  
  <pre><code>size = 10000
  bins = 30
  intuitive_models = {"model1":(2, 5),
                      "model2":(1, 4, 1),
                      "model3":(0.8, 4, 1, 1.25),
                      "fitted":truncnorm.fit(values)}
  
  # store the tuncnorm random sample into a dict
  model_results = dict()
  for model, params in intuitive_models.items():
      model_results[model] = truncnorm(*params).rvs(size)
  
  # plot the random sample vs the oserved values
  for model, params in intuitive_models.iteritems():
      plt.figure()
      plt.hist(model_results[model], bins=bins, normed=True)
      plt.title(model + ': ' + repr(params))
      plt.hist(values, normed=True, alpha=0.5)
  
  # tabular comparison    
  print pd.DataFrame(model_results).describe()
  </code></pre>
  
  <p>which produced the following tabular data:</p>
  
  <pre><code>             fitted        model1        model2        model3 
  count  10000.000000  10000.000000  10000.000000  10000.000000 
  mean       1.024707      2.372819      2.524923      2.698601 
  std        0.014362      0.333144      0.443857      0.584215 
  min        1.000019      2.000040      2.000007      2.000019 
  25%        1.012248      2.121838      2.181642      2.245088 
  50%        1.024518      2.280975      2.407814      2.557983 
  75%        1.036996      2.534782      2.757778      2.998948 
  max        1.049991      4.829619      4.982337      5.905201
  </code></pre>
  
  <p>Thanks Bertie.
  p.s. I hope this is accepted as a coding question and not a stats question - which is why I have posted it here.</p>
  
  <p><img src="https://i.stack.imgur.com/Ropff.png" alt="model1">
  <img src="https://i.stack.imgur.com/mz1S8.png" alt="model2">
  <img src="https://i.stack.imgur.com/yUw1P.png" alt="model3">
  <img src="https://i.stack.imgur.com/LfBzu.png" alt="Fitted"></p>
  
  <p>-- Update 28-Aug-2014 --
  The idea with this post was to see get some help with the scipy.stats.truncnorm.fit method and in the couple of days, I have built my own clunky algorithm. From my discussions with Robert, I get the impression the R or standard implementation for truncnorm is to only take 3 parameters. For those coming to this post later, once scipy has an improved fitting engine, this is what I have estimated (assuming we want an asymptotic right tail). </p>
  
  <p><img src="https://i.stack.imgur.com/r299Y.png" alt="enter image description here"></p>