How can I impute a large multivariate dataset? I am researching a medical dataset that has a large number and complexity of variables (> 900 variabes, > 3000 with dummy variables). Not only biomarkers and the patient's demographical data is included, but also questionnaires (continuous, non-ordered categorical and ordered categorical, non-respectively). Due to this complexity I have not managed to use MICE in R or Python, it errors out and states message that the data is numerically impossible to solve.
This lead me to try to impute the questionnaires apart from the other variables, because they are the largest part of the dataset and because the correlation between missing data has the highest probability. Because it is a longitudinal study, some people might not have been able to fill out the first questionnaire for example, because they were not yet included in the study, so by using the other questionnaires we might be able to impute for the missing data of the first questionnaire. 
The method I tried to use to impute for the questionnaires was MICE, which failed, again. The other method I tried was order the questions according to missing data, where v1 had the least amount of missing data and vn the largest amount. I then tried to find variables that had no overlapping missing data with v1, let's call them [ v2 ... vn]. I then tried to create a logistic regression between variables [v2 ... vn] after deleting the patients with missing data, the problem was that the logistic regression could not converge, because only 1 of the two classes was available (there were no 1's available in the data after deleting the patients with missing data for the predictor variables: [v2 ... vn])
Because this also failed I am currently going to try to do this: find the 10% most correlating questions with v1 that do not have any overlapping data with v1: [v2 ... vn], impute these predictors ([v2 ... vn]) with their mean and then logistically regress them with the non-missing data in v1. Then I will predict the missing data in v1 with the non-missing data from [v2 ... vn]. Now, even though this is likely to work, it does not seem very statistically correct and that is my main question, is this a solid method to use for imputing a large multivariate dataset, if not, what can I improve or are there other methods that I have overlooked?
 A: The algorithm you have described at the bottom is basically one "chaining" iteration of MICE [1], so it is less crazy than you'd guess...
From your description, you have a lot of columns and they are of mixed type (categorical, numeric). Now, if the proportion of missing values is large per column, then it is close to hopeless to let an algorithm solve that problem. On the other hand, if the amout of missing values is typically small per column (5% or so), then you can give a try to MICE with e.g. random forest as model technique. Why? Because it can deal with categorical (not dummy coded!) and numeric fiels. Furthermore, they work quite well without tedious hyperparameter tuning.
Then, depending on the purpose of the analysis, you can use multiple imputed versions of the data set, run your analysis on each of them and then pool the results of the analyses (this is the tricky part of multiple imputation). 
Let me illustrate with my own R version [2] of imputation by random forests, a fast alternative to [3] based on the very efficient ranger package [4]:
# install.packages(devtools)
library(devtools) 
install_github("mayer79/missRanger")

library(missRanger)

# Generate data with missing values in all columns
set.seed(39)
head(irisWithNA <- generateNA(iris))

# Show data set with missings
      Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1          5.1         3.5          1.4         0.2  setosa
2          4.9         3.0          1.4         0.2  setosa
3          4.7         3.2          1.3         0.2  setosa
4          4.6          NA          1.5         0.2  setosa
5          5.0         3.6          1.4         0.2    <NA>
6          5.4         3.9          1.7          NA  setosa

# Impute missing values 
head(irisImputed <- missRanger(irisWithNA, pmm.k = 3, seed = 75757))

# Gives
  Sepal.Length Sepal.Width Petal.Length Petal.Width Species
1          5.1         3.5          1.4         0.2  setosa
2          4.9         3.0          1.4         0.2  setosa
3          4.7         3.2          1.3         0.2  setosa
4          4.6         3.7          1.5         0.2  setosa
5          5.0         3.6          1.4         0.2  setosa
6          5.4         3.9          1.7         0.4  setosa

Just run with different seeds to get realisticly varying versions of the same data. The "predictive mean matching" aspect is crucial if you are doing multiple imputation. 
[1] Van Buuren, S., Groothuis-Oudshoorn, K. (2011). mice: Multivariate Imputation by Chained Equations in R. Journal of Statistical Software, 45(3), 1-67. http://www.jstatsoft.org/v45/i03/
[2] https://github.com/mayer79/missRanger
[3] Stekhoven, D.J. and Buehlmann, P. (2012), 'MissForest - nonparametric missing value imputation for mixed-type data', Bioinformatics, 28(1) 2012, 112-118, doi: 10.1093/bioinformatics/btr597
[4] Wright, M. N. & Ziegler, A. (2016). ranger: A Fast Implementation of Random Forests for High Dimensional Data in C++ and R. Journal of Statistical Software, in press. http://arxiv.org/abs/1508.04409.
