def expected_information_gain(self, expparams):
        r"""
        Calculates the expected information gain for each hypothetical experiment.

        :param expparams: The experiments at which to compute expected
            information gain.
        :type expparams: :class:`~numpy.ndarray` of dtype given by the current
            model's :attr:`~qinfer.abstract_model.Simulatable.expparams_dtype` property,
            and of shape ``(n,)``

        :return float: The expected information gain for each 
            hypothetical experiment in ``expparams``.
        """
        # This is a special case of the KL divergence estimator (see below),
        # in which the other distribution is guaranteed to share support.
        
        # for models whose outcome number changes with experiment, we 
        # take the easy way out and for-loop over experiments
        n_eps = expparams.size
        if n_eps > 1 and not self.model.is_n_outcomes_constant:
            risk = np.empty(n_eps)
            for idx in range(n_eps):
                risk[idx] = self.expected_information_gain(expparams[idx, np.newaxis])
            return risk
        
        # number of outcomes for the first experiment
        os = self.model.domain(expparams[0,np.newaxis])[0].values

        # compute the hypothetical weights, likelihoods and normalizations for
        # every possible outcome and expparam
        # the likelihood over outcomes should sum to 1, so don't compute for last outcome
        w_hyp, L, N = self.hypothetical_update(
                os[:-1], 
                expparams, 
                return_normalization=True, 
                return_likelihood=True
            )
        w_hyp_last_outcome = (1 - L.sum(axis=0)) * self.particle_weights[np.newaxis, :]
        N = np.concatenate([N[:,:,0], np.sum(w_hyp_last_outcome[np.newaxis,:,:], axis=2)], axis=0)
        w_hyp_last_outcome = w_hyp_last_outcome / N[-1,:,np.newaxis]
        w_hyp = np.concatenate([w_hyp, w_hyp_last_outcome[np.newaxis,:,:]], axis=0)
        # w_hyp.shape == (n_out, n_eps, n_particles)
        # N.shape == (n_out, n_eps)

        # compute the Kullback-Liebler divergence for every experiment and possible outcome
        # KLD.shape == (n_out, n_eps)
        KLD = np.sum(w_hyp * np.log(w_hyp / self.particle_weights), axis=2)

        # return the expected KLD (ie expected info gain) for every experiment
        return np.sum(N * KLD, axis=0)
r"""
        Calculates the expected information gain for each hypothetical experiment.

        :param expparams: The experiments at which to compute expected
            information gain.
        :type expparams: :class:`~numpy.ndarray` of dtype given by the current
            model's :attr:`~qinfer.abstract_model.Simulatable.expparams_dtype` property,
            and of shape ``(n,)``

        :return float: The expected information gain for each 
            hypothetical experiment in ``expparams``.