C:/programs/etirm/src/swig_etirm.cpp

Go to the documentation of this file.

/*! \file swig_etirm.cpp
 
  \brief
  Definitions of functions and classes to generate scripting language wrapper 
  functions for ETIRM using SWIG.

  The file swig_etirm_types.h (which is included in swig_etirm.h) must be 
  provided which should contain typedef's for the following types:
 
  ResponseVector  Vector containing item responses.
  item_type   Item class.
  examinee_type Examinee class.
  lvdist_type   Latent variable distribution class.
  random_type   Type of random number generator used for bootstrap.

  For example:

  // Do not process with SWIG, SWIG cannot handle the template arguments
  xx_ifndef SWIG
 
  // Vector of item responses
  typedef std::vector<etirm::Response> ResponseVector;
 
  // Class to hold information about discrete latent variable distribution
  typedef etirm::DiscreteLatentDist<etirm::Real> lvdist_type;
 
  // Class to hold information about examinees
  typedef etirm::ExamineeGrpCov<ResponseVector, etirm::RealVector> examinee_type;
 
  // Type used for array of items. Using ItemNR allows different items to be modeled
  // by different classes descendent from ItemNR (e.g., the 3PL model could be used for
  // some items, and the 2PL model used for other items).
  typedef etirm::ItemNR<lvdist_type> item_type;
 
  // Use the Mersenne Twister (http://www.math.keio.ac.jp/matumoto/emt.html)
  // from the Boost random number library
  // (http://www.boost.org/libs/random/index.html) 
  // as random number generator for bootstrap samples.
  typedef boost::mt19937 random_type;

  xx_endif // SWIG
 
  An application using these wrapper functions should declare a subclass
  of SwigEtirmRun whose constructor initializes the item object pointers, and 
  assigns them to the 'items' data member of SwigEtirmRun. 
  An application should provide functions for creating and deleting an object
  that is a subclass of SwigEtirmRun which is stored in the global variable
  gEtirmRun defined in this file.

  Definitions of the function CheckRunInit must be provided. This
  function is declared in swig_etirm.h, but each application must define it.
  It is not defined in swig_etirm.cpp.
 
  The function CheckRunInit checks whether gEtirmRun has been
  initialized, and if not throws an exception. The message contained
  in the exception will differ for different applications. 

  Estimation Toolkit for Item Response Models (ETIRM)
  http://www.smallwaters.com/software/cpp/etirm.html

  Author(s): 
  Werner Wothke, maintenance (http://www.smallwaters.com)
  Brad Hanson (http://www.b-a-h.com/)
  See the file LICENSE for information on usage and redistribution.

  Copyright (C) 2008, Werner Wothke
  Copyright (c) 2000-2002, Bradley A. Hanson
 */

#ifdef ETIRM_NO_DIR_PREFIX
#include "swig_etirm.h"
#include "MStepIRT.h"
#include "DiscreteNormalDist.h"
#include "ItemParamPriorBeta4.h"
#include "ItemParamPriorLogNormal.h"
#include "ItemParamPriorNormal.h"
#include "BootstrapSample.h"
#include "SimulateResponses.h"
#include "ExamineeThetaMLE.h"
#else
#include "etirm/swig_etirm.h"
#include "etirm/MStepIRT.h"
#include "etirm/DiscreteNormalDist.h"
#include "etirm/ItemParamPriorBeta4.h"
#include "etirm/ItemParamPriorLogNormal.h"
#include "etirm/ItemParamPriorNormal.h"
#include "etirm/BootstrapSample.h"
#include "etirm/SimulateResponses.h"
#include "etirm/ExamineeThetaMLE.h"
#endif

// Use integer random number generator from boost library
// (http://www.boost.org).
#include <boost/random/uniform_int.hpp>

#include <cstdio>  // sprintf prototype
// for compilers which do not put C library functions in std namespace
#ifdef BOOST_NO_STDC_NAMESPACE
namespace std
{ using ::sprintf;}
#endif

// Definitions not wrapped by SWIG
namespace etirm
{

  //! Global variable holding information about a run.
  SwigEtirmRun *gEtirmRun = 0;

  /*!
    \brief
    Increments counts based on item response.
   
    \section function_args Function Parameters
   
    \param[in]  respIndex "One" offset index corresponding to response to increment 
        count for (1 = first response category, 2 = second response category, etc.).
    \param[in]  count   Amount to increment counts.
    \param[in]  group   Group to increment counts for (1, 2, ...).
   */
  void ItemRespCounts::AddResponse(int respIndex, Real count, int group)
  {

    if (respIndex > catCounts.size() || respIndex < 1)
    {
      throw InvalidArgument("Invalid response index", "ItemRespCounts::AddResponse");
    }

    totalCount[0] += count; // total count
    totalCount[group] += count; // count for group

    catCounts(1, respIndex) += count; // total count
    catCounts(group+1, respIndex) += count; // count for group
  }

  /*!
    \brief
    Class constructor.
   
    \section function_args Function Parameters
   
    \param[in]  nitems    Total number of items.
    \param[in]  nlatentcat  Number of categories of discrete latent variable.
    \param[in]  ngroups   Number of examinee groups.
    \param[in]  minTheta  Minimum value of theta points for discrete latent 
        variable distribution.
    \param[in]  maxTheta  Maximum value of theta points for discrete latent 
        variable distribution.
    \param[in]  uniquePoints  If true then unique latent distribution points 
        are used for each examinee group.
   */
  SwigEtirmRun::SwigEtirmRun(int nitems, int nlatentcat, int ngroups, Real minTheta, Real maxTheta,
      bool uniquePoints) :
    numItems(nitems), numLatentCat(nlatentcat), numGroups(ngroups), examineeCounts(ngroups+1, 0.0),
        items(nitems), itemStats(nitems), minProc(nitems), latentDist(nlatentcat, ngroups,
            uniquePoints), rand_boot(0), base_rand_simulate(0), rand_simulate(0)
  {
    int i;

    // Set default points and weights for latent variable distribution to standard normal
    // for first group
    DiscreteNormalDist(nlatentcat, minTheta, maxTheta, latentDist.begin_points(1),
        latentDist.begin_weights(1));

    // Set points and weights for other groups equal to the points and weights for group 1
    for (i=2; i<=ngroups; ++i)
    {
      DiscreteLatentDist<Real>::weight_iterator w1 = latentDist.begin_weights(1);
      DiscreteLatentDist<Real>::weight_iterator w2 = latentDist.begin_weights(i);
      for (int j = nlatentcat; j--; ++w1, ++w2)
      {
        *w2 = *w1;
      }

      if (uniquePoints)
      {
        // initialize unique points for examinee group
        DiscreteLatentDist<Real>::point_iterator w1p = latentDist.begin_points(1);
        DiscreteLatentDist<Real>::point_iterator w2p = latentDist.begin_points(i);
        for (int j = nlatentcat; j--; ++w1p, ++w2p)
        {
          *w2p = *w1p;
        }
      }
    }

    /* initialize vectors of pointers to null */
    for (i=0; i<nitems; ++i)
    {
      items[i] = 0;
      itemStats[i] = 0;
      minProc[i] = 0;
    }

  }

  //! Releases memory allocated by constructor.
  SwigEtirmRun::~SwigEtirmRun()
  {
    int i;

    for (i=0; i<numItems; i++)
    {
      if (items[i])
      {
        items[i]->DeletePriors();
        delete items[i];
      }

      if (itemStats[i])
        delete itemStats[i];
      if (minProc[i])
        delete minProc[i];
    }

    int n = examinees.size();
    for (i=0; i<n; ++i)
    {
      delete examinees[i];
    }

    delete rand_boot;
    delete base_rand_simulate;
    delete rand_simulate;
  }

  // Check that data for examinees exist
  void SwigEtirmRun::CheckExaminees(const char *funcname)
  {
    if ((gEtirmRun->examinees).size() == 0)
    {
      throw RuntimeError("Error: No examinee data", funcname);
    }
  }

  // Check that examinee number is valid
  void SwigEtirmRun::CheckExamineeNumber(int examno, const char *funcname)
  {
    if ( (examno < 1) || ( examno > (gEtirmRun->examinees).size() ) )
    {
      char errstr[50];
      std::sprintf(errstr, "Invalid examinee number: %d", examno);
      gEtirmRun->returnString = errstr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), funcname);
    }
  }

  // Check that item number is valid
  void SwigEtirmRun::CheckItemNumber(int itemno, const char *funcname)
  {
    if (itemno < 1 || itemno >(gEtirmRun->numItems) ) // Syntax edited, ww, 2-24-2008.
    {
      char errstr[50];
      std::sprintf(errstr, "Invalid item number: %d", itemno);
      gEtirmRun->returnString = errstr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), funcname);
    }
  }

  // Check that group number is valid
  //  group   Value to check.
  //  funcname  Name of calling function (used in error message)
  void SwigEtirmRun::CheckGroup(int group, const char *funcname)
  {
    if (group < 1 || group >(gEtirmRun->numGroups) ) // Syntax edited, ww, 2-24-2008
    {
      char errstr[50];
      std::sprintf(errstr, "Invalid examinee group: %d", group);
      gEtirmRun->returnString = errstr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), funcname);
    }
  }

  /*! 
    \brief
    Tests whether item parameter index is valid for the item.
   
    \section function_args Function Parameters
    
    \param[in]  *item   Pointer to item object for which index is checked.
    \param[in]  index   Zero-offset index of item parameter.
    \param[in]  *funcname pointer to name of calling function (used in error message).
   */
  void CheckItemParam(item_type *item, int index, const char *funcname)
  {
    if (item->NumParameters() <= index || index < 0)
    {
      char errstr[100];
      std::sprintf(errstr, "Invalid item parameter index for item %d: %d", item->Index()+1, index);
      gEtirmRun->returnString = errstr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), funcname);
    }
  }

  /*!
    \brief
    Create a prior distribution object and return a pointer to it.
   
    \section function_args Function Parameters

    \param[in]  pstr  Type of prior distribution, valued as ("none", "normal", "lognormal", "beta").
    \param[in]  priorparam  Parameters of prior distribution.
    \param[in]  funcname  Name of calling function (used in error messages)
   */
  ItemParamPrior *CreatePrior(const std::string &pstr, const double_vector &priorparam,
      const char *funcname)
  {
    const char *err_message = "Invalid number of prior parameters";

    if (pstr.compare("none") == 0)
    {
      return 0;
    }
    else if (pstr.compare("normal") == 0)
    {
      if (priorparam.size() != 2)
        throw RuntimeError(err_message, funcname);
      return new ItemParamPriorNormal(priorparam[0], priorparam[1]);
    }
    else if (pstr.compare("lognormal") == 0)
    {
      if (priorparam.size() != 2)
        throw RuntimeError(err_message, funcname);
      return new ItemParamPriorLogNormal(priorparam[0], priorparam[1]);
    }
    else if (pstr.compare("beta") == 0)
    {
      if (priorparam.size() != 4)
        throw RuntimeError(err_message, funcname);
      return new ItemParamPriorBeta4(priorparam[0], priorparam[1], priorparam[2], priorparam[3]);
    }
    else
    {
      throw RuntimeError("Invalid prior type", funcname);
    }
    return 0;
  }

  /*!
    \brief
    Converts a response for an item into a character, where '0' represents the 
    first response, '1' represents the second response, etc.
   
    If the number of response categories for the item is greater than 10, then the 
    characters returned will be ascii characters greater than '9'. For example, a 
    ':' is return for a response in response category 11, since ':' is one greater
    than '9' in the ascii sequence.

    \section function_args Function Parameters
    
    \param[in]  r Item response to be converted.
    \param[in]  *item Pointer to item object.
   */
  char Resp2Char(Response r, const item_type *item)
  {
    const char zero = '0';
    Response np = item_type::NotPresentedResponse();
    char cr;
    if (r == np)
      cr = np - item->FirstResponse() + zero;
    else
      cr = item->ResponseIndex(r) + zero;

    return cr;
  }

  /****** Functions for which SWIG wrappers are generated ******/

  /* Member functions for estep class */

  /*!
    \brief
    Initialize estep object.
   
    \section function_args Function Parameters
   
    \param[in]  *itemno Pointer to list of item numbers to use for computing 
    examinee posteriors distribution in E-Step. If itemno = NUL, use all items (default).
   */
  estep::estep(int_vector *itemno)
  {
    CheckRunInit();

    if (itemno)
    {
      mItems = new ItemVector(ItemSubset(itemno,gEtirmRun->items,"new_estep"));

      mEStep = new estep_type(mItems->begin(), mItems->end(), gEtirmRun->latentDist);
    }
    else
    {
      mItems = 0;
      mEStep = new estep_type((gEtirmRun->items).begin(), (gEtirmRun->items).end(), gEtirmRun->latentDist);
    }
  }

  /*!
    \brief
    Release memory allocated in constructor.
   */
  estep::~estep()
  {
    if (mItems)
      delete mItems;
    delete mEStep;
  }

  /*!
    \brief
    Calls mEStep->DoEStep to perform E-step.
   
    Returns marginal loglikelihood of examinees' responses (sum over examinees of the 
    marginal loglikelihood of an examinee's responses) plus sum of prior likelihoods 
    over all item parameters.  This is the value of the marginal posterior density that 
    the EM algorithm is maximizing at the values of the item parameters computed in the 
    last M-step.
   
    \section function_args Function Parameters
   
    \param[in]  compute_post  Flag: If compute_post == TRUE, then compute posterior 
    distributions for all examinees. If compute_post == FALSE, then use previously 
    stored posteriors for examinees. (Default: TRUE).
   
    \param[in]  store_post  Flag: If store_post == TRUE, then store the posterior 
    distribution computed for each examinee. These posterior distributions are 
    stored as part of the examinee objects.
    (Default: FALSE).
   
    \param[in]  *estep_items  Pointer to list of items for which n and r are updated. If a 
    null pointer is passed then n arnd r are updated for all items used in the E-step 
    to compute examinee posterior distributions.
   */
  double estep::compute(bool compute_post, bool store_post, int_vector *estep_items) // compute_post and store_post retyped as "bool", ww, 2-24-2008.
  {
    const char *fname = "estep_compute";
    CheckRunInit();
    gEtirmRun->CheckExaminees(fname);

    if (estep_items)
    {
      if (estep_items->size() > 0)
      {
        ItemVector itemsub = ItemSubset(estep_items, gEtirmRun->items, fname);
        return mEStep->DoEStep((gEtirmRun->examinees).begin(), (gEtirmRun->examinees).end(), itemsub.begin(), itemsub.end(), compute_post, store_post);
      }
      else // do not update n and r for any items
      {
        return mEStep->DoEStep((gEtirmRun->examinees).begin(), (gEtirmRun->examinees).end(), (gEtirmRun->items).begin(), (gEtirmRun->items).begin(), compute_post, store_post);
      }
    }
    else
    {
      return mEStep->DoEStep((gEtirmRun->examinees).begin(), (gEtirmRun->examinees).end(), compute_post, store_post);
    }
  }

  /*!
    \brief
    Assigns missing response code to identify examinees who did not respond to an item.
   */
  void set_missing_resp(char nr)
  {
    if (gEtirmRun)
    {
      throw RuntimeError("The not presented response can only be set before any items are initialized", 0);
    }
    Item<Real>::SetNotPresentedResponse(nr);
  }

  /*!
    \brief
    Returns the number of items.
   */
  int num_items()
  {
    CheckRunInit();

    return gEtirmRun->numItems;
  }

  /*!
    \brief
    Returns number of categories of the discrete theta distribution.
   */
  int num_latent_dist_points()
  {
    CheckRunInit();

    return gEtirmRun->numLatentCat;
  }

  /*!
    \brief
    Returns number of examinee groups.
   */
  int num_groups()
  {
    CheckRunInit();

    return gEtirmRun->numGroups;
  }

  /*!
    \brief
    Returns number of examinees.
   */
  int num_examinees()
  {
    CheckRunInit();

    return (gEtirmRun->examinees).size();
  }

  /*!
    \brief
    Returns the name of model used for an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno Item number (1-offset).
   */
  const char * item_get_model(int itemno)
  {

    CheckRunInit();
    std::string &modelName = gEtirmRun->returnString;
    gEtirmRun->CheckItemNumber(itemno, "item_get_model");

    item_type *item = (gEtirmRun->items)[itemno-1];
    modelName = item->ModelName();

    return modelName.c_str();
  }

  /*!
    \brief
    Assigns a value to one item parameter of an item.
   
    \section function_args Function Parameters
   
    \param[in] paramno  1-offset index of parameter in parameter vector for item.
    \param[in] itemno   Item number (1-offset).
    \param[in] paramvalue Value parameter is set to.
   */
  void item_set_param(int paramno, int itemno, double paramvalue)
  {
    const char *fname = "item_set_param";
    CheckRunInit();
    int index = paramno - 1;

    gEtirmRun->CheckItemNumber(itemno, fname);
    item_type *item = (gEtirmRun->items)[itemno-1];
    CheckItemParam(item, index, fname);

    /* Check that value of parameter has nonzero prior density */
    item_type::prior_iterator pri = item->PriorsIterator() + index;
    if (*pri && (*pri)->ZeroDensity(paramvalue))
    {
      char merr[100];
      std::sprintf(merr, "Parameter specified for item %d has zero prior density", itemno);
      gEtirmRun->returnString = merr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), fname);
    }

    item_type::param_iterator pv = item->ParametersIterator();

    pv[index] = paramvalue;
  }

  /*!
    \brief
    Assigns values to all item parameters of an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
    \param[in]  *params  Pointer to params vector containing the values to set the 
        parameters to.
   
    Note: The order of the parameters in the params vector is: 
    a (if there is an a parameter), one or more b's, and
    c (if there is a c parameter).
   */
  void item_set_params(int itemno, double_vector *params)
  {
    const char *fname = "item_set_params";
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, fname);

    item_type *item = (gEtirmRun->items)[itemno-1];

    // Assign new parameters
    item_type::param_iterator pv = item->ParametersIterator();
    double_vector::iterator iv = params->begin();
    int n;
    for (n = item->NumParameters(); n--; ++pv, ++iv)
    {
      *pv = *iv;
    }

    /* Check that values of all parameters have nonzero prior density */
    iv = params->begin();
    item_type::prior_iterator pri = item->PriorsIterator();
    for (n = item->NumParameters(); n--; ++pri, ++iv)
    {
      if (*pri && (*pri)->ZeroDensity(*iv))
      {
        char merr[100];
        std::sprintf(merr, "Parameter specified for item %d has zero prior density", itemno);
        throw InvalidArgument(merr, fname);
      }
    }
  }

  /*!
    \brief
    Assigns values to all fixed and estimated item parameters of an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
    \param[in]  *params  Pointer to params vector containing the values to set the 
        parameters to.
   
    Note: The order of the parameters in the params vector is: 
    a (if there is an a parameter), one or more b's, and
    c (if there is a c parameter).
   */
  void item_set_all_params(int itemno, double_vector *params)
  {
    const char *fname = "item_set_all_params";
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, fname);

    item_type *item = (gEtirmRun->items)[itemno-1];

    // Assign new fixed and estimated parameters
    item->SetAllParameters(*params);

    /* Check that values of estimated parameters all have nonzero prior density */
    item_type::param_iterator iv = item->ParametersIterator();
    item_type::prior_iterator pri = item->PriorsIterator();
    for (int n = item->NumParameters(); n--; ++pri, ++iv)
    {
      if (*pri && (*pri)->ZeroDensity(*iv))
      {
        char merr[100];
        std::sprintf(merr, "Parameter specified for item %d has zero prior density", itemno);
        gEtirmRun->returnString = merr;
        throw InvalidArgument((gEtirmRun->returnString).c_str(), fname);
      }
    }
  }

  /*!
    \brief
    Returns value of one item parameter of an item.
   
    \section function_args Function Parameters
    
    \param[in]  paramno  1-offset index of the parameter in the item's parameter vector.
    \param[in]  itemno  Item number (1-offset).
   
    Note: The order of the parameters in the params vector is: 
    a (if there is an a parameter), one or more b's, and
    c (if there is a c parameter).
   */
  double item_get_param(int paramno, int itemno)
  {
    const char *fname = "item_get_param";
    CheckRunInit();
    int index = paramno - 1;

    gEtirmRun->CheckItemNumber(itemno, fname);
    item_type *item = (gEtirmRun->items)[itemno-1];
    CheckItemParam(item, index, fname);

    item_type::param_iterator pv = item->ParametersIterator();

    return pv[index];
  }

  /*!
    \brief
    Returns vector of all estimated item parameters values of an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
   
    Note: The order of the parameters in the params vector is: 
    a (if there is an a parameter), one or more b's, and
    c (if there is a c parameter).
   */
  double_vector *item_get_params(int itemno)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "item_get_params");

    item_type *item = (gEtirmRun->items)[itemno-1];

    int n = item->NumParameters();
    double_vector *params = new double_vector(n);

    item_type::param_iterator pv = item->ParametersIterator();
    double_vector::iterator iv = params->begin();
    for (; n--; ++pv, ++iv)
    {
      *iv = *pv;
    }

    return params;
  }

  /*!
    \brief
    Returns vector of all fixed and estimated item parameters values of an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
   
    Note: The order of the parameters in the params vector is: 
    a (if there is an a parameter), one or more b's, and
    c (if there is a c parameter).
   */
  double_vector *item_get_all_params(int itemno)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "item_get_all_params");

    item_type *item = (gEtirmRun->items)[itemno-1];

    RealVector allparam = item->GetAllParameters();
    int n = allparam.size();
    double_vector *params = new double_vector(n);

    RealVector::iterator pv = allparam.begin();
    double_vector::iterator iv = params->begin();
    for (; n--; ++pv, ++iv)
    {
      *iv = *pv;
    }

    return params;
  }

  /*!
    \brief
    Returns the number of parameters of an item.

    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
   */
  int item_num_params(int itemno)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "item_num_params");

    item_type *item = (gEtirmRun->items)[itemno-1];

    return item->NumParameters();
  }

  /*!
    \brief
    Returns the number of response categories of an item.

    \section function_args Function Parameters
   
    \param[in]  itemno  Item number (1-offset).
   */
  int item_num_resp_cat(int itemno)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "item_num_resp_cat");

    item_type *item = (gEtirmRun->items)[itemno-1];

    return item->NumRespCat();
  }

  /*!
    \brief
    Assigns prior distribution parameters for one item parameter.

    \section function_args Function Parameters
   
    \param[in]  paramno 1-offset index of parameter in parameter vector for item.
    \param[in]  itemno  Number of item for which prior is set (1-offset).
    \param[in]  *priortype Pointer to type of prior distribution ("normal", "lognormal", "beta", "none").
    \param[in]  *dlist Pointer to vector of prior distribution parameters.
   */
  void item_set_prior(int paramno, int itemno, char *priortype, double_vector *dlist)
  {
    const char *funcname = "item_set_prior";
    CheckRunInit();
    int index = paramno - 1;

    gEtirmRun->CheckItemNumber(itemno, funcname);
    item_type *item = (gEtirmRun->items)[itemno-1];
    CheckItemParam(item, index, funcname);

    item_type::prior_iterator pv = item->PriorsIterator() + index;

    if (*pv)
      delete *pv;

    std::string pstr(priortype);
    *pv = CreatePrior(pstr, *dlist, funcname);
  }

  /*!
    \brief
    Returns type of prior distribution ("normal", "lognormal", "beta", "none")
    for one item parameter.

    \section function_args Function Parameters
   
    \param[in]  paramno 1-offset index of parameter in parameter vector for item.
    \param[in]  itemno  Number of item for which prior is returned (1-offset).
   */
  const char *item_get_prior_type(int paramno, int itemno)
  {
    const char *fname = "item_get_prior_type";
    CheckRunInit();
    std::string &priorName = gEtirmRun->returnString;
    int index = paramno - 1;

    gEtirmRun->CheckItemNumber(itemno, fname);
    item_type *item = (gEtirmRun->items)[itemno-1];
    CheckItemParam(item, index, fname);

    item_type::prior_iterator pv = item->PriorsIterator() + index;

    if (*pv != 0)
    {
      priorName = (*pv)->DistributionName();
    }
    else
    {
      priorName = "none";
    }

    return priorName.c_str();
  }

  /*!
    \brief
    Returns vector of prior distribution parameters for one item parameter.

    \section function_args Function Parameters
   
    \param[in]  paramno 1-offset index of parameter in parameter vector for item.
    \param[in]  itemno  Number of item for which prior is returned (1-offset).
   */
  double_vector *item_get_prior_param(int paramno, int itemno)
  {
    const char *fname = "item_get_prior_param";
    CheckRunInit();
    int index = paramno - 1;

    gEtirmRun->CheckItemNumber(itemno, fname);
    item_type *item = (gEtirmRun->items)[itemno-1];
    CheckItemParam(item, index, fname);

    item_type::prior_iterator pv = item->PriorsIterator() + index;

    double_vector *v;
    if (*pv != 0)
    {
      int n = (*pv)->NumParameters();

      v = new double_vector(n);

      *v = (*pv)->GetParameters();
    }
    else
    {
      v = new double_vector();
    }

    return v;
  }

  /*!
    \brief
    Returns vector of response counts in each response category of an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Number of item for which response counts are returned (1-offset).
    \param[in]  group Group to return counts for (1-offset). Specify group = 0 to return 
        counts across all groups.
   */
  double_vector *item_cat_counts(int itemno, int group)
  {
    const char *fname = "item_cat_counts";
    CheckRunInit();

    gEtirmRun->CheckItemNumber(itemno, fname);
    if (group != 0)
      gEtirmRun->CheckGroup(group, "item_resp_count");

    ItemRespCounts *stats = (gEtirmRun->itemStats)[itemno-1];

    double_vector *v = new double_vector(stats->CategoryCounts(group));

    return v;
  }

  /*!
    \brief
    Returns the number of examinees responding to an item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Number of item for which response counts are returned (1-offset).
    \param[in]  group Group to return counts for (1-offset). Specify group = 0 to return 
        count across all groups.
   */
  double item_resp_count(int itemno, int group)
  {
    CheckRunInit();

    gEtirmRun->CheckItemNumber(itemno, "item_resp_count");
    if (group != 0)
      gEtirmRun->CheckGroup(group, "item_resp_count");

    return (gEtirmRun->itemStats)[itemno-1]->RespCount(group);
  }
  /*!
    \brief
    Transforms the parameter estimates of an item to a different latent variable scale.
   
    Returns zero if parameters are successfully scaled, or returns nonzero if
    scaling would result in an invalid parameter value.   

    \section function_args Function Parameters
   
    \param[in]  itemno  Number of item for which to transform parameters (1-offset).
    \param[in]  slope  Slope of scale transformation to apply to item parameters.
    \param[in]  intercept Intercept of scale transformation to apply to item parameters.
    \param[in]  ignorePriorError  If true then do not report an error if transformed parameter
        has zero density in prior used for that parameter.
   */
  int item_scale_params(int itemno, double slope, double intercept, bool ignorePriorError) 
  { // Retyped ignorePriorErro from "int" to "bool", ww, 2-24-2008.
    CheckRunInit();

    gEtirmRun->CheckItemNumber(itemno, "item_scale_params");
    item_type *item = (gEtirmRun->items)[itemno-1];

    return item->ScaleParameters(slope, intercept, ignorePriorError);
  }

  /*!
    \brief
    Returns the probability that an examinee with a particular theta value will give a 
    particular response to a particular item.
   
    \section function_args Function Parameters
   
    \param[in]  itemno  Number (1-based) of item for which response probability is returned.
    \param[in]  response  Integer representing response, where a response in the first
        response category is 0, a response in the second response category is 1, etc.
    \param[in]  theta Value of latent variable for which response probability is returned.
   */
  double item_prob_resp(int itemno, int response, double theta)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "item_prob_resp");

    item_type *item = (gEtirmRun->items)[itemno-1];

    return item->ProbResp(item->IndexResponse(response), theta);
  }

  /*!
    \brief
    Returns a vector of values of the test characteristic curve, given a vector of theta values and
    a selection of items. 
   
    It is assumed for each item that the lowest response category corresponds to a score of 0, the
    second response category corresponds to a score of 1, etc.

    \section function_args Function Parameters
   
    \param[in]  *thetas Pointer to vector of theta values.
    \param[in]  *ilist2 Pointer to vector of item sequence numbers (1-based).
   */
  double_vector *test_characteristic_curve(double_vector *thetas, int_vector *ilist2)
  {
    CheckRunInit();

    double_vector tcc(thetas->size(), 0.0);

    // Vector to hold response probabilities for an item.
    double_vector probs;

    ItemVector *sitems;
    if (ilist2)
    {
      sitems = new ItemVector(ItemSubset(ilist2, gEtirmRun->items, "test_characteristic_curve"));
    }
    else
    {
      sitems = &(gEtirmRun->items);
    }

    // loop over items
    ItemVector::iterator ii = sitems->begin();
    for (int i = sitems->size(); i--; ++ii)
    {
      int ncat = (*ii)->NumRespCat();
      probs.newsize(ncat);

      double_vector::iterator icc = tcc.begin();
      double_vector::iterator it = thetas->begin();
      for (int j = tcc.size(); j--; ++icc, ++it) // loop over thetas
      {
        (*ii)->ProbRespAll(*it, probs.begin());

        // first response category corresponds to 0, second response category
        // corresponds to 1, etc.
        double dresp = 1.0; // Start with response 1, 0 is skipped
        double sum = 0.0;
        double_vector::iterator ip = probs.begin()+1;
        for (int k = ncat-1; k--; ++ip, ++dresp)
        {
          sum += dresp * *ip;
        }

        // Add expected score for item to TCC
        *icc += sum;
      }
    }

    if (ilist2)
      delete sitems;

    return new double_vector(tcc);
  }

  /*!
    \brief
    Assigns quadrature point values to each discrete category of the discrete latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  *dlist Pointer to vector of quadrature points of the latent variable distribution.
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  void dist_set_points(double_vector *dlist, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_set_points");

    int n = dlist->size();
    if (n != (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Number of points does not match number of categories in latent distribution",
          "dist_set_points");
    }

    DiscreteLatentDist<Real>::point_iterator ip = (gEtirmRun->latentDist).begin_points(group);

    double_vector::iterator id = dlist->begin();
    for (; n--; ++ip, ++id)
    {
      *ip = *id;
    }
  }

  /*!
    \brief
    Assigns a quadrature point value to one discrete category of the discrete latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  index Index of latent variable category (1, 2, ..., number of categories).
    \param[in]  p Value to assign to point.
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  void dist_set_point(int index, double p, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_set_point");

    if (index < 1 || index> (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Invalid category of latent variable distribution","dist_set_point");
    }

    DiscreteLatentDist<Real>::point_iterator ip = (gEtirmRun->latentDist).begin_points(group);

    ip[index-1] = p;
  }

  /*!
    \brief
    Returns vector of quadrature point values of the discrete latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  double_vector *dist_get_points(int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_get_points");

    DiscreteLatentDist<Real>::point_iterator ip = (gEtirmRun->latentDist).begin_points(group);

    int n = (gEtirmRun->latentDist).size();
    double_vector *outv = new double_vector(n);
    double_vector::iterator id = outv->begin();
    for (; n--; ++ip, ++id)
    {
      *id = *ip;
    }

    return outv;
  }

  /*!
    \brief
    Returns the quadrature point value of one discrete category of the discrete latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  index Index of latent variable category (1, 2, ..., number of categories).
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  double dist_get_point(int index, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_get_point");

    if (index < 1 || index> (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Invalid category of latent variable distribution","dist_get_point");
    }

    DiscreteLatentDist<Real>::point_iterator ip = (gEtirmRun->latentDist).begin_points(group);

    return ip[index-1];
  }

  /*!
    \brief
    Assigns quadrature weights to the quadrature points of the latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  *dlist Pointer to vector of quadrature weights of the latent variable distribution.
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  void dist_set_probs(double_vector *dlist, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_set_probs");

    int n = dlist->size();
    if (n != (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Number of probabilities does not match number of categories in latent distribution",
          "dist_set_probs");
    }

    // assign probabilities
    DiscreteLatentDist<Real>::weight_iterator ip = (gEtirmRun->latentDist).begin_weights(group);
    double_vector::iterator id = dlist->begin();
    Real sum = 0.0;
    for (; n--; ++ip, ++id)
    {
      *ip = *id;
      sum += *ip;
    }

    // standardize probabilities to they sum to 1
    ip = (gEtirmRun->latentDist).begin_weights(group);
    for (n = dlist->size(); n--; ++ip)
    {
      *ip /= sum;
    }
  }

  /*!
    \brief
    Assigns the quadrature weight to one discrete category of the discrete latent variable distribution.

    \section function_args Function Parameters
   
    \param[in]  index Index of latent variable category (1, 2, ..., number of categories).
    \param[in]  w Value to assign to weight.
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  void dist_set_prob(int index, double w, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_set_prob");

    if (index < 1 || index> (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Invalid category of latent variable distribution","dist_set_prob");
    }

    DiscreteLatentDist<Real>::weight_iterator ip = (gEtirmRun->latentDist).begin_weights(group);

    ip[index-1] = w;
  }

  /*!
    \brief
    Returns vector of quadrature weights of the latent variable distribution for one group of
    examinees.

    \section function_args Function Parameters
   
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  double_vector *dist_get_probs(int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_get_probs");

    DiscreteLatentDist<Real>::weight_iterator ip = (gEtirmRun->latentDist).begin_weights(group);

    int n = (gEtirmRun->latentDist).size();
    double_vector *outv = new double_vector(n);
    double_vector::iterator id = outv->begin();
    for (; n--; ++ip, ++id)
    {
      *id = *ip;
    }

    return outv;

  }

  /*!
    \brief
    Returns the quadrature weight of one discrete category of the latent variable 
    distribution for one group of examinees.

    \section function_args Function Parameters
   
    \param[in]  index Index of latent variable category (1, 2, ..., number of categories).
    \param[in]  group Group to use (1, 2, ..., number of groups), default = 1.
   */
  double dist_get_prob(int index, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_get_prob");

    if (index < 1 || index> (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Invalid category of latent variable distribution","dist_get_prob");
    }

    DiscreteLatentDist<Real>::weight_iterator ip = (gEtirmRun->latentDist).begin_weights(group);

    return ip[index-1];
  }

  /*!
    \brief
    Transforms the quadrature points of latent variable distribution to a new scale.

    \section function_args Function Parameters
   
    \param[in]  slope Slope parameter of linear scale transformation.
    \param[in]  intercept Intercept parameter of linear scale transformation.
   */
  void dist_transform(double slope, double intercept)
  {
    CheckRunInit();

    (gEtirmRun->latentDist).Transform(slope, intercept);
  }

  /*!
    \brief
    Scales to the quadrature points of latent variable distribution to yield a 
    specfic mean and standard deviation in one group.

    \section function_args Function Parameters
   
    \param[in]  mean  Mean of reference group after scaling the latent distribution.
    \param[in]  sd  Standard deviation of reference group after scaling the latent distribution.
    \param[in]  group Group to use as reference group (1, 2, ..., number of groups).
   */
  double_vector *dist_scale(double mean, double sd, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_scale");

    double slope, intercept;
    (gEtirmRun->latentDist).Scale(mean, sd, group, slope, intercept);

    double_vector *outv = new double_vector(2);

    (*outv)[0] = slope;
    (*outv)[1] = intercept;

    return outv;
  }

  /*!
    \brief
    Returns a vector with the mean and standard deviation of the latent variable distribution 
    for a selected group of examinees.

    \section function_args Function Parameters
   
    \param[in]  group Group to use as reference group (1, 2, ..., number of groups).
   */
  double_vector *dist_mean_sd(int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "dist_mean_sd");

    double mean, sd;
    (gEtirmRun->latentDist).MeanSD(group, mean, sd);

    double_vector *outv = new double_vector(2);

    (*outv)[0] = mean;
    (*outv)[1] = sd;

    return outv;
  }

  /*!
    \brief
    Returns 1 if unique sets of quadrature points are used for two or more groups of 
    examinees, or returns 0 otherwise.
   */
  int dist_unique_points()
  {
    CheckRunInit();

    return ((gEtirmRun->latentDist).NumGroupsUnique() > 1 && (gEtirmRun->latentDist).NumGroups() > 1);
  }

  /*!
    \brief
    Returns a vector of probabilities for a discrete distribution to approximate a normal distribution 
    over a set of equally spaced points.

    Note: The mean and s.d. parameters only effect the quadrature points themselves, so it does not
    matter which values are passed for the mean and s.d.

    \section function_args Function Parameters
   
    \param[in]  npoints   Number of discrete points in the latent distribution.
    \param[in]  minPoint  Minimum point of the discrete distribution (must be smaller than the mean).
    \param[in]  maxPoint  Maximum point of the discrete distribution (must be larger than the mean).
    \param[in]  mean  Mean of the normal distribution.
    \param[in]  sd  Standard deviation of the normal distribution.
   */
  double_vector *normal_dist_prob(int npoints, double minPoint, double maxPoint, double mean,
      double sd)
  {
    RealVector points(npoints);
    RealVector *prob = new RealVector(npoints);

    DiscreteNormalDist(npoints, minPoint, maxPoint, points.begin(), prob->begin(), mean, sd);

    return prob;
  }

  /*!
    \brief
    Returns a vector of quadrature points for a discrete distribution over a set of equally-spaced 
    points.
    \section function_args Function Parameters
   
    \param[in]  npoints   Number of discrete points in the latent distribution.
    \param[in]  minPoint  Minimum point of the discrete distribution (must be smaller than the mean).
    \param[in]  maxPoint  Maximum point of the discrete distribution (must be larger than the mean).
    \param[in]  mean  Mean of the normal distribution.
    \param[in]  sd  Standard deviation of the normal distribution.
   */
  double_vector *normal_dist_points(int npoints, double minPoint, double maxPoint, double mean,
      double sd)
  {
    RealVector prob(npoints);
    RealVector *points = new RealVector(npoints);

    // Probabilities do not matter so make up mean and sd
    DiscreteNormalDist(npoints, minPoint, maxPoint, points->begin(), prob.begin(), mean, sd);

    return points;
  }
  
  /*!
    \brief
    CalculateS M-step for a set of items.

    \section function_args Function Parameters
   
    \param[in]  ignore_max_iter Flag: If ignore_max_iter == TRUE, then exceeding the maximum number of 
        iterations in optimization procedure is NOT a fatal error; if ignore_max_iter == FALSE, then the
        exceeding the max_iter is treated as a fatal error (default = FALSE).
    \param[in]  *itemno  Pointer to vector of integers giving item numbers (1-offset)
       for which the M-step is calculated. If itemno is the null pointer, then all items are used
       (default = null).
   */
  int mstep_items(bool ignore_max_iter, int_vector *itemno)
  {  // Retyped ignore_max_iter from "int" to "bool", ww, 2-25-2008.
    CheckRunInit();

    if (itemno)
    {
      const char *fname = "mstep_items";
      ItemVector itemsub = ItemSubset(itemno, gEtirmRun->items, fname);
      UncminVector minsub = ItemSubset(itemno, gEtirmRun->minProc, fname);
      return MStepItems(minsub.begin(), itemsub.begin(), itemsub.end(), gEtirmRun->mstepMaxDiff,
          ignore_max_iter);
    }
    else
    {
      return MStepItems((gEtirmRun->minProc).begin(), (gEtirmRun->items).begin(), (gEtirmRun->items).end(), gEtirmRun->mstepMaxDiff, ignore_max_iter);
    }
  }

  /*!
    \brief
    Returns message from M-step minimization for one item.

    \section function_args Function Parameters
   
    \param[in]  itemno  Integer giving item number (1-offset).
   */
  int mstep_message(int itemno)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "mstep_message");

    return ((gEtirmRun->minProc)[itemno-1])->GetMessage();
  }

  /*!
    \brief
    Assigns the maximum number of iterations permitted with optimization procedure
    in mstep_items for one item.
      
    \section function_args Function Parameters
   
    \param[in]  itemno  Integer giving item number (1-offset).
    \param[in]  maxiter  Maximum number of iterations for stepwise optimization.
   */
  void mstep_max_iter(int itemno, int maxiter)
  {
    CheckRunInit();
    gEtirmRun->CheckItemNumber(itemno, "mstep_max_iter");

    return ((gEtirmRun->minProc)[itemno-1])->SetMaxIter(maxiter);
  }

  /*!
    \brief
    Returns the maximum relative difference between parameter estimates in two
    successive EM iterations, as computed in last call to mstep_items.
   */
  double mstep_max_diff()
  {
    CheckRunInit();

    return gEtirmRun->mstepMaxDiff;
  }

  /*!
    \brief
    Executes M-step for the latent distribution in one group.
    
    Returns maximum relative difference between new and old probabilities.
      
    \section function_args Function Parameters
   
    \param[in]  *e  Pointer to not sure what(?!). To-Do: Find out what *e relates to! 
    \param[in]  group Number of examinee group (1, 2, ..., number of groups).
   */
  double mstep_dist(estep *e, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckGroup(group, "mstep_dist");

    estep_type::ngroup_iterator i = (e->GetEStep())->GetNGroup(group);
    RealVector::size_type n = (e->GetEStep())->size();
    RealVector prob(i, i+n);
    return (gEtirmRun->latentDist).MStep(prob, group);
  }

  /*!
    \brief
    Add an examinee objec to the end of the examinee vector.
    
    Returns the examinee number corresponding to this examinee.

    \section function_args Function Parameters
  
    \param[in]  *responses  Pointer to vector of integer-formatted item responses. 
        A negative integer indicates the examinee did not respond to the item.
        The integers representing responses are zero-offset (A response in the 
        first response category is represented as zero).
    \param[in]  group   Integer from 1 to maximum number of examinee groups
        giving the group the examinee belongs to.
    \param[in]  count   Frequency given to response pattern (usually, count = 1).
   */
  int add_examinee(int_vector *responses, int group, double count)
  {
    const char *fname = "add_examinee";
    CheckRunInit();
    gEtirmRun->CheckGroup(group, fname);

    int len = responses->size();

    if (len != gEtirmRun->numItems)
    {
      throw RuntimeError("Invalid number of item responses", fname);
    }

    examinee_type *e = new examinee_type(len, group);

    ResponseVector r(len);
    ItemVector::iterator iitem = (gEtirmRun->items).begin();
    int_vector::iterator iir = responses->begin();
    ResponseVector::iterator ir = r.begin();
    ItemStatsVector::iterator is = (gEtirmRun->itemStats).begin();
    Response np = Item<Real>::NotPresentedResponse();
    for (; len--; ++iitem, ++iir, ++ir, ++is)
    {
      if (*iir < 0)
        *ir = np;
      else
      {
        *ir = (*iitem)->IndexResponse(*iir);
        if (!((*iitem)->ValidResponse(*ir)))
        {
          throw RuntimeError("Invalid item response", fname);
        }
        (*is)->AddResponse(*iir+1, count, group);
      }
    }

    e->SetResponses(r);
    e->SetCount(count);

    (gEtirmRun->examinees).push_back(e);

    (gEtirmRun->examineeCounts)[0] += count; // total count
    (gEtirmRun->examineeCounts)[group] += count; // group count

    return (gEtirmRun->examinees).size();
  }

  /*!
    \brief
    Returns pointer to an integer vector of examinee item responses to all items.

    A non-negative integer indicates an examinee response in the category 
    represented by the integer (0 = first response category, 1 = second response 
    category, etc.). A negative integer indicates the examinee did not respond to 
    the item.
   
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee for whom item responses are returned.
   */
  int_vector *examinee_responses(int examineeno)
  {
    const char *fname = "examinee_responses";
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, fname);

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    if ((gEtirmRun->items).size() != examinee->NumItems())
    {
      throw RuntimeError("Number of examinee responses does not match number of items", fname);
    }

    // iterators to first and one past last response
    examinee_type::response_iterator first = examinee->responses_begin();
    examinee_type::response_iterator last = examinee->responses_end();

    int_vector *responses = new int_vector(examinee->NumItems());
    int_vector::iterator ir = responses->begin();
    ItemVector::iterator ii = (gEtirmRun->items).begin();
    Response np = Item<Real>::NotPresentedResponse();
    while (first != last)
    {
      if (*first == np)
        *ir = -1;
      else
        *ir = (*ii)->ResponseIndex(*first);
      ++ir;
      ++ii;
      ++first;
    }

    return responses;
  }

  /*!
    \brief
    Returns string containing (character) examinee responses to all items.
    
    A response in the first response category is '0', a response in the second category is 
    '1', etc. Because each response is represented by a single character care must be 
    taken when some items have more than 10 response categories. For example, a response 
    in the 11th response category would be represented by the ascii character ':', which 
    is the 10th character greater than '0' in the ascii character sequence.
   
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee for whom item responses are returned.
   */

  const char *examinee_response_str(int examineeno)
  {
    const char *fname = "examinee_response_str";
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, fname);

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    if ((gEtirmRun->items).size() != examinee->NumItems())
    {
      throw RuntimeError("Number of examinee responses does not match number of items", fname);
    }

    // iterators to first and one past last response
    examinee_type::response_iterator first = examinee->responses_begin();
    examinee_type::response_iterator last = examinee->responses_end();

    std::string &responseStr = gEtirmRun->returnString;
    responseStr.clear();
    ItemVector::iterator ii = (gEtirmRun->items).begin();
    while (first != last)
    {
      responseStr.push_back(Resp2Char(*first, *ii));
      ++ii;
      ++first;
    }

    return responseStr.c_str();
  }

  /*!
    \brief
    Return the number of the group the examinee belongs to.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee for whom group membership is to be returned.
   */

  int examinee_get_group(int examineeno)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_group");

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    return examinee->Group();
  }

  /*!
    \brief
    Assigns which group the examinee belongs to.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee for whom group membership is to be assigned.
    \param[in]  group   Group number to be assigned to examinee.
   */
  void examinee_set_group(int examineeno, int group)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_set_group");

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    return examinee->SetGroup(group);
  }

  /*!
    \brief
    Assigns the count (or weight) of the examinee.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose count (or weight) is to be assigned.
    \param[in]  count   Count to be assigned to examinee.
   */
  void examinee_set_count(int examineeno, double count)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_set_count");

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    examinee->SetCount(count);
  }

  /*!
    \brief
    Returns the count (or weight) of the examinee.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose count (or weight) is to be returned.
   */
  double examinee_get_count(int examineeno)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_get_count");

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    return examinee->Count();
  }

  /*!
    \brief
    Assigns posterior distribution of an examinee.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose posterior is to be assigned.
    \param[in]  *posterior  Pointer to vector of posterior probabilities of the examinee.
   */
  void examinee_set_posterior(int examineeno, double_vector *posterior)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_set_posterior");

    if (posterior->size() != (gEtirmRun->latentDist).size())
    {
      throw RuntimeError("Invalid number of posterior probabilities",
          "examinee_set_posterior");
    }

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];

    examinee_type::posterior_vector epost(gEtirmRun->numLatentCat);
    examinee_type::posterior_vector::iterator iep = epost.begin();
    double_vector::iterator ip = posterior->begin();

    // assign probabilities
    int i;
    Real sum = 0.0;
    for (i = gEtirmRun->numLatentCat; i--; ++iep, ++ip)
    {
      *iep = *ip;
      sum += *iep;
    }

    // standardize probabilities to sum to 1.0
    iep = epost.begin();
    for (i = gEtirmRun->numLatentCat; i--; ++iep)
    {
      *iep /= sum;
    }

    examinee->SetPosterior(epost);
  }

  /*!
    \brief
    Returns pointer to vector of the posterior distribution of an examinee.
    
    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose posterior is to be returned.
   */ 
  double_vector* examinee_get_posterior(int examineeno)
  {
    char *fname = "examinee_get_posterior";
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, fname);

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];
    int n = examinee->NumLatentVarCat();

    if (n == 0)
    {
      char merr[100];
      std::sprintf(merr, "Posterior has not been computed for examinee %d", examineeno);
      gEtirmRun->returnString = merr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), fname);
    }

    double_vector *post = new double_vector(n);

    examinee_type::posterior_vector::iterator ie = examinee->posterior_begin();
    double_vector::iterator ip = post->begin();
    for (; n--; ++ie, ++ip)
    {
      *ip = *ie;
    }

    return post;
  }

  /*!
    \brief
    Returns means of an examinee's posterior distribution.
    
    Posterior must have already been computed, for example, by estep_compute.

    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose posterior mean is to be returned.
   */ 
  double examinee_posterior_mean(int examineeno)
  {
    char *fname = "examinee_posterior_mean";
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, fname);

    examinee_type *examinee = (gEtirmRun->examinees)[examineeno-1];
    int n = examinee->NumLatentVarCat();

    if (n == 0)
    {
      char merr[100];
      std::sprintf(merr, "Posterior distribution has not been computed for examinee %d", examineeno);
      gEtirmRun->returnString = merr;
      throw RuntimeError((gEtirmRun->returnString).c_str(), fname);
    }

    examinee_type::posterior_vector::iterator iw = examinee->posterior_begin();
    lvdist_type::point_iterator ip = (gEtirmRun->latentDist).begin_points();
    double mean = 0.0;
    for (; n--; ++iw, ++ip)
    {
      mean += *ip * *iw;
    }

    return mean;
  }

  /*!
    \brief
    Returns an examinee's maximum likelihood estimate of theta.
    
    Posterior must have already been computed, for example, by estep_compute.

    \section function_args Function Parameters
  
    \param[in]  examineeno  Number of examinee whose theta MLE is to be returned (1 = first examinee).
    \param[in]  minTheta  Minimum value of theta estimate.
    \param[in]  maxTheta  Maximum value of theta estimate.  
    \param[in]  precision Length of interval in which MLE is determined to lie. This should be greater 
        than, roughly, 3.0e-8.
    \param[in]  *itemno5    Pointer to a list of item numbers of the items to use in computing examinee MLE. 
        If itemno5 = NULL, then all items will be used.
   */
  double examinee_theta_MLE(int examineeno, double minTheta, double maxTheta, double precision,
      int_vector *itemno5)
  {
    CheckRunInit();
    gEtirmRun->CheckExamineeNumber(examineeno, "examinee_theta_MLE");

    if (!gEtirmRun->base_rand_simulate)
    {
      gEtirmRun->base_rand_simulate = new random_type();
      gEtirmRun->rand_simulate = new boost::uniform_01<random_type>(*(gEtirmRun->base_rand_simulate));
    }

    ItemVector *sitems;
    if (itemno5)
    {
      sitems = new ItemVector(ItemSubset(itemno5, gEtirmRun->items, "examinee_theta_MLE"));
    }
    else
    {
      sitems = &(gEtirmRun->items);
    }

    double theta =
        ExamineeThetaMLE<ItemVector::iterator, examinee_type::response_vector::iterator>(minTheta,
            maxTheta, precision, sitems->begin(), sitems->end(), (gEtirmRun->examinees[examineeno-1])->responses_begin());

    if (itemno5)
      delete sitems;

    return theta;
  }

  /*!
    \brief
    Returns total examinee count in an examinee group (1, 2, ...), or across all groups.
    
    \section function_args Function Parameters
  
    \param[in]  group Group to return count for (1-offset). Specify group = 0 to return 
        counts across all groups.
   */
  double examinees_count(int group)
  {
    CheckRunInit();
    if (group != 0)
      gEtirmRun->CheckGroup(group, "examinee_count");

    return (gEtirmRun->examineeCounts)[group];
  }

  /*!
    \brief
    Assigns seed of random number generator used for bootstrap samples.
    
    \section function_args Function Parameters
  
    \param[in]  seed  Seed value, 0 <= seed <= 4,294,967,295.
   */
  void bootstrap_seed(unsigned long seed)
  {
    CheckRunInit();

    // convert to type used for seed in boost random number library
    boost::uint32_t boost_seed = seed;

    if (gEtirmRun->rand_boot)
    {
      gEtirmRun->rand_boot->seed(boost_seed);
    }
    else
    {
      gEtirmRun->rand_boot = new random_type(boost_seed);
    }
  }

  /*!
    \brief
    Generate bootstrap sample of examinees.
   */
  void bootstrap_sample()
  {
    CheckRunInit();
    gEtirmRun->CheckExaminees("bootstrap_sample");

    if (!gEtirmRun->rand_boot)
      gEtirmRun->rand_boot = new random_type();

    boost::uniform_int<random_type> urand(*(gEtirmRun->rand_boot), 1, (gEtirmRun->examinees).size());

    BootstrapSample((gEtirmRun->examinees).begin(), (gEtirmRun->examinees).end(), urand);

    // Set vector of examinee counts
    gEtirmRun->examineeCounts = 0.0; // initialize to zero
    ExamineeVector::iterator ie = (gEtirmRun->examinees).begin();
    for (int i = (gEtirmRun->examinees).size(); i--; ++ie)
    {
      double count = (*ie)->Count();
      (gEtirmRun->examineeCounts)[0] += count; // total count
      int group = (*ie)->Group();
      (gEtirmRun->examineeCounts)[group] += count; // group count
    }
  }

  /*!
    \brief
    Assigns seed of random number generator used for simulating item responses.
    
    \section function_args Function Parameters
  
    \param[in]  seed  Seed value, 0 <= seed <= 4,294,967,295.
   */  
  void simulate_seed(unsigned long seed)
  {
    CheckRunInit();

    // convert to type used for seed in boost random number library
    boost::uint32_t boost_seed = seed;

    if (gEtirmRun->base_rand_simulate)
    {
      gEtirmRun->base_rand_simulate->seed(boost_seed);
    }
    else
    {
      gEtirmRun->base_rand_simulate = new random_type(boost_seed);
      gEtirmRun->rand_simulate = new boost::uniform_01<random_type>(*(gEtirmRun->base_rand_simulate));
    }
  }

  /*!
    \brief
    Simulates item responses (integers) for a specific value of the latent variable.  

    Returns pointer to integer vector containing item responses, where
    item responses are integers from 0 to one minus the maximum
    number of response categories for the item.
    
    \section function_args Function Parameters
  
    \param[in]  theta Value of latent variable for which item responses are generated.
    \param[in]  itemno2 List of item numbers of the items for which responses will be
        simulated.
   */
  int_vector *simulate_responses(double theta, int_vector *itemno2)
  {
    CheckRunInit();

    if (!gEtirmRun->base_rand_simulate)
    {
      gEtirmRun->base_rand_simulate = new random_type();
      gEtirmRun->rand_simulate = new boost::uniform_01<random_type>(*(gEtirmRun->base_rand_simulate));
    }

    int_vector *resp;
    if (itemno2)
    {
      resp = new int_vector(itemno2->size());
      ItemVector sitems = ItemSubset(itemno2, gEtirmRun->items, "simulate_responses");
      resp->newsize(sitems.size());
      SimulateResponses(sitems.begin(), sitems.end(), theta, *(gEtirmRun->rand_simulate),
          resp->begin(), true);
    }
    else
    {
      resp = new int_vector((gEtirmRun->items).size());
      SimulateResponses((gEtirmRun->items).begin(), (gEtirmRun->items).end(), theta, *(gEtirmRun->rand_simulate), resp->begin(), true);
    }

    return resp;
  }

  /*!
    \brief
    Simulates item responses (characters) for a specific value of the latent variable.  

    Returns simulated responses in a string, where each character of the string gives a 
    response to one item. The first response category of an item is represented by a 
    '0', the second response category by a '1', etc. 

    Because each response is represented by a single character care must be taken when 
    some items have more than 10 response categories. For example, a response in the 
    11th response category would be represented by the ascii character ':', which is 
    the 10th character greater than '0' in the ascii character sequence.
       
    \section function_args Function Parameters
  
    \param[in]  theta Value of latent variable for which item responses are generated.
    \param[in]  itemno2 List of item numbers of the items for which responses will be
        simulated.
   */
  const char *simulate_response_str(double theta, int_vector *itemno2)
  {
    CheckRunInit();

    // Vector of integers representing simulated responses
    int_vector *iresp = simulate_responses(theta, itemno2);

    // Resize string to contain item responses
    std::string &responseStr = gEtirmRun->returnString;
    responseStr.resize(iresp->size());

    // Convert integers to characters
    int_vector::iterator ii = iresp->begin();
    std::string::iterator ic = responseStr.begin();
    for (int i = iresp->size(); i--; ++ii, ++ic)
    {
      if (*ii < 0)
      {
        *ic = item_type::NotPresentedResponse();
      }
      else
      {
        *ic = *ii + '0';
      }
    }

    delete iresp;

    return responseStr.c_str();
  }

  /*!
    \brief
    Reads item responses from an input record.  

    The item responses are assumed to be integers from 0 to one minus the maximum
    number of response categories for the item. A response representing a missing 
    response is assigned -1.

    Returns poiner to an integer vector containing thew item responses.
       
    \section function_args Function Parameters
  
    \param[in]  *line  Pointer to character string to read the item responses from.
    \param[in]  *offset  Pointer to vector of zero-based offsets of each item response in 'line'.
    \param[in]  *len   Pointer to vector of field widths for each item response.  
   */
  int_vector *get_responses(char *line, int_vector *offset, int_vector *len)
  {
    const char *fname = "get_responses";
    CheckRunInit();

    int n = offset->size();

    if (n != len->size())
    {
      throw InvalidArgument("Lengths of offset and length vectors do not match", fname);
    }

    int_vector *resp = new int_vector(n);

    int_vector::iterator ir = resp->begin();
    int_vector::iterator ioff = offset->begin();
    int_vector::iterator ilen = len->begin();
    Response np = Item<Real>::NotPresentedResponse();
    for (; n--; ++ioff, ++ilen, ++ir)
    {
      // Convert response to integer
      char *pos = line + *ioff + *ilen - 1;
      *ir = 0;
      int power = 1;
      for (int i=*ilen; i--; power *= 10, --pos)
      {
        if (*pos == np) // check for missing item response
        {
          *ir = -1;
          break;
        }
        else if ( (*pos < '0') && (*pos > '9')) // check that item response is valid // Sytax edited, ww, 2-24-2008.
        {
          char errstr[50];
          std::sprintf(errstr, "Invalid item response: %1c", *pos);
          gEtirmRun->returnString = errstr;
          throw RuntimeError((gEtirmRun->returnString).c_str(), fname);
        }
        *ir += power * (*pos - '0');
      }
    }

    return resp;
  }

  /*!
    \brief
    Read item responses from a string where responses to only some items are
    present. The responses to the remaining items are assumed to be missing.

    The item responses are assumed to be integers from 0 to one minus the maximum
    number of response categories for the item. A response representing a missing 
    response is assigned -1.

    Returns poiner to an integer vector containing thew item responses.
       
    \section function_args Function Parameters
  
    \param[in]  *line  Pointer to character string to read the item responses from.
    \param[in]  *offset  Pointer to vector of zero-based offsets of each item response in 'line'.
    \param[in]  *len   Pointer to vector of field widths for each item response.  
    \param[in]  *items  Pointer to vector 1-offset indices of items for which responses are read.
   */
  int_vector *get_responses_missing(char *line, int_vector *offset, int_vector *len,
      int_vector *items)
  {
    const char *fname = "get_responses_missing";
    CheckRunInit();

    int n = offset->size();

    if (n != len->size() || n != items->size())
    {
      throw InvalidArgument("Lengths of vector arguments do not match", fname);
    }

    int_vector *resp = get_responses(line, offset, len);

    // initially assign all missing responses
    int nall = gEtirmRun->numItems;
    int_vector *allresp = new int_vector(nall, -1);

    // Assign responses read
    int_vector::iterator ir = resp->begin();
    int_vector::iterator ii = items->begin();
    for (; n--; ++ii, ++ir)
    {
      if ( (*ii < 1) || (*ii > nall)) // Syntax edited, ww, 2-24-2008.
      {
        delete resp;
        throw InvalidArgument("Invalid item number", fname);
      }
      (*allresp)[*ii-1] = *ir;
    }

    delete resp;

    return allresp;
  }

} // namespace etirm

---