Implement `LogisticPCA` as a distinct variant of matrix decomposition useful for binary data

### Describe the workflow you want to enable

Currently, there is no included implementation of a PCA algorithm made for handling binary data in the scikit-learn library. However, the algorithm for "logistic PCA" is well founded, although different methods for its estimation exist, and it meets the inclusion criteria of being at least 3 years since first publication, having over 200 citations between known papers, and has a wide use and usefulness. The proposed implementation follows that of Lee et al. (2010), which allows for the use of an optional L1 regularization term. The intent is to create a `LogisticPCA` implementation that mimics the existing `SparsePCA` implementation to the greatest degree possible acknowledging their important differences, and which utilizes the adopted best practices to fit within the existing matrix decomposition algorithms and the current API of scikit-learn.

References

Tipping, Michael E. "Probabilistic visualisation of high-dimensional binary data." Advances in Neural Information Processing Systems (1999): 592-598.

Lee, Seokho, Jianhua Z. Huang, and Jianhua Hu. "Sparse logistic principal components analysis for binary data." The Annals of Applied Statistics, 4.3 (2010): 1579-1601.

Landgraf, Andrew J. and Yoonkyung Lee. "Dimensionality reduction for binary data through the projection of natural parameters." Journal of Multivariate Analysis, v.180 (2020): 104668.

### Describe your proposed solution

The sigmoid function, defined as `sig = lambda x: 1 / (1 + np.exp(-x))`, plays a crucial role in the logistic PCA algorithm. Its purpose is to map any real-valued number into the range (0, 1), which is particularly useful for binary data. 

Here are the roles the sigmoid function plays in the algorithm:

Probability Estimation: The sigmoid function is used to estimate probabilities. In the context of logistic PCA, it helps in modeling the probability that a given binary variable is 1. This is essential for handling binary data, where each element in the data matrix `X` is either 0 or 1.

Variational Parameters: The lambda function `lam = lambda x: (0.5 - sig(x)) / (2 * x)` uses the sigmoid function to compute a term that is part of the variational parameters. These parameters are crucial for approximating the posterior distribution of the latent variables.

Gradient and Hessian Computation: During the optimization steps, the sigmoid function is used to compute gradients and Hessians. These are necessary for updating the model parameters (`W` and `b`) to minimize the reconstruction error.

By transforming the data through the sigmoid function, the logistic PCA algorithm can effectively handle the binary nature of the input data, making it suitable for applications where the data is represented as binary variables.

Here are the key steps involved:

Initialization:

- Initialize the identity matrix `I`, the projection matrix `W`, latent variables `mu`, bias `b`, covariance matrices `C`, and variational parameters `xi`. Using the constants `n_components`, `n_samples`, and `n_features` where applicable in their definitions.

- Define the sigmoid function `sig` and the lambda function `lam`.

Iterative Optimization:

The optimization process runs for a specified number of iterations to `max_iter`.

Step 1. Update Sufficient Statistics:
- For each sample `n`:
   - Extract the sample `x_n` from the data matrix `X`.
   - Compute the approximation `lam_n` using the lambda function `lam`.
   - Update the covariance matrix `C` and the latent variables `mu` using the computed values.

Step 2. Optimize Variational Parameters:
- For each sample `n`:
   - Compute the posterior statistics `z` and `E_zz`.
   - Calculate `xixi`, which involves the projection matrix `W`, the posterior statistics, and the bias `b`.
   - Update the variational parameters `xi` based on the computed values.

Step 3. Update Model Parameters:
- Compute the expected value `E_zhzh` for each sample `n`.
- Compute the expected value `E_zh`.
- For each dimension `i`:
   - Compute the approximation `lam_i` using the lambda function `lam`.
   - Calculate the gradient `g` and the Hessian `H`. Regularization can be applied with the optional `lasso_lambda` parameter.
   - Solve the linear system to update the projection matrix `W` and the bias `b`.

Return Results:
- After completing the specified number of iterations, return the estimated projection matrix `W`, latent variables `mu`, and bias `b`.

### Describe alternatives you've considered, if relevant

_No response_

### Additional context

In another project of mine, I have a working implementation of logistic PCA based on Mikael Brudfors's (2021) implementation under an MIT License, which I had modified to include an L1 regularization parameter per the recommendation of Lee et al. (2010). However, it is my intention, and I've already started that process, to rework the code for greater efficiency and clarity. Given that, the above-described solution will be altered as necessary while maintaining the same essential functionality. Also, I want to make sure that I don't cut Mikael out of this process, especially since to my knowledge that code was the first known implementation of this algorithm in Python. I'm seeking their involvement if they desire to be a part of this addition to the scikit-learn ecosystem if it is tentatively approved for inclusion as a new feature.

COMMENT:
@GaelVaroquaux Personally, I've used logistic PCA as an embedding technique following the one-hot encoding of categorical data. Frequently, as many other people likely encounter in everyday situations, I'm wanting to use mixed categorical and numerical datasets. Outside of something like Gower's distance, I've had the need at times to use other distance metrics robust for such data to create functional coordinate projections inside of the developing topological data analysis and topological machine learning frameworks. 

But, as I suggested above, another method I've employed are turning to various embedding techniques to collapse those "broken out" binary features back down into something that a wider variety of distance metrics can be computed on. Obviously, this is an area where a lot of new work is currently being done, and it may not be the typical use for something like this. Mainly, I'm doing my own exploratory investigations of comparing and contrasting between things like Cat2Vec and logistic PCA for use along with these topological methods.

However, in doing that, it came to mind that having at least one known implementation of a logistic matrix decomposition with optional regularization for inducing sparsity within the scikit-learn ecosystem would be helpful to others as well for different purposes. As PCA is a commonly used dimensionality reduction technique, it struck me as a good candidate to suggest.

If this is too niche of an ask right now, I understand if it might be best for me (and anyone else who might be interested in helping along the way) to move forward with it as an interoperable external module.

Labels: New Feature, Needs Decision - Include Feature