Beberapa imputasi untuk nilai yang hilang

Saya ingin menggunakan imputasi untuk mengganti nilai yang hilang dalam kumpulan data saya di bawah batasan tertentu.

Sebagai contoh, saya ingin variabel yang diperhitungkan x1menjadi lebih besar atau sama dengan jumlah dari dua variabel saya yang lain, katakan x2dan x3. Saya juga ingin x3diperhitungkan oleh salah satu 0atau >= 14dan saya ingin x2diperhitungkan oleh salah satu 0atau >= 16.

Saya mencoba mendefinisikan batasan ini dalam SPSS untuk beberapa imputasi, tetapi dalam SPSS saya hanya dapat menentukan nilai maksimum dan minimum. Apakah ada cara untuk mendefinisikan kendala lebih lanjut dalam SPSS atau Anda tahu paket R yang akan membiarkan saya mendefinisikan kendala tersebut untuk penghilangan nilai yang hilang?

Data saya adalah sebagai berikut:

   x1 =c(21, 50, 31, 15, 36, 82, 14, 14, 19, 18, 16, 36, 583, NA,NA,NA, 50, 52, 26, 24)
   x2 = c(0, NA, 18,0, 19, 0, NA, 0, 0, 0, 0, 0, 0,NA,NA, NA, 22, NA, 0, 0)
   x3 = c(0, 0, 0, 0, 0, 54, 0 ,0, 0, 0, 0, 0, 0, NA, NA, NA, NA, 0, 0, 0)
   dat=data.frame(x1=x1, x2=x2, x3=x3)
   > dat
       x1 x2 x3
   1   21  0  0
   2   50 NA  0
   3   31 18  0
   4   15  0  0
   5   36 19  0
   6   82  0 54
   7   14 NA  0
   8   14  0  0
   9   19  0  0
   10  18  0  0
   11  16  0  0
   12  36  0  0
   13 583  0  0
   14  NA NA NA
   15  NA NA NA
   16  NA NA NA
   17  50 22 NA
   18  52 NA  0
   19  26  0  0
   20  24  0  0

Salah satu solusinya adalah menulis fungsi imputasi kustom Anda sendiri untuk micepaket tersebut. Paket ini disiapkan untuk ini dan pengaturannya tanpa rasa sakit.

Pertama-tama kita mengatur data seperti yang disarankan:

dat=data.frame(x1=c(21, 50, 31, 15, 36, 82, 14, 14, 19, 18, 16, 36, 583, NA,NA,NA, 50, 52, 26, 24), 
               x2=c(0, NA, 18,0, 19, 0, NA, 0, 0, 0, 0, 0, 0,NA,NA, NA, 22, NA, 0, 0), 
               x3=c(0, 0, 0, 0, 0, 54, 0 ,0, 0, 0, 0, 0, 0, NA, NA, NA, NA, 0, 0, 0))

Selanjutnya kita memuat micepaket dan melihat metode apa yang dipilih secara default:

library(mice)
# Do a non-imputation
imp_base <- mice(dat, m=0, maxit = 0)

# Find the methods that mice chooses
imp_base$method
# Returns: "pmm" "pmm" "pmm"

# Look at the imputation matrix
imp_base$predictorMatrix
# Returns:
#   x1 x2 x3
#x1  0  1  1
#x2  1  0  1
#x3  1  1  0

The pmmsingkatan rata pencocokan prediktif - mungkin algoritma imputasi paling populer untuk imputing variabel kontinu. Ini menghitung nilai prediksi menggunakan model regresi dan memilih 5 elemen terdekat dengan nilai prediksi (dengan jarak Euclidean ). Elemen-elemen yang dipilih disebut kolam donor dan nilai akhir dipilih secara acak dari kolam donor ini.

Dari matriks prediksi kami menemukan bahwa metode mendapatkan variabel lulus yang menarik bagi pembatasan. Perhatikan bahwa baris adalah variabel target dan kolom prediktornya. Jika x1 tidak memiliki 1 di kolom x3 kita harus menambahkan ini dalam matriks:imp_base$predictorMatrix["x1","x3"] <- 1

Sekarang ke bagian yang menyenangkan, menghasilkan metode imputasi. Saya telah memilih metode yang agak kasar di sini di mana saya membuang semua nilai jika tidak memenuhi kriteria. Ini dapat mengakibatkan waktu loop yang panjang dan mungkin berpotensi lebih efisien untuk menjaga imputasi yang valid dan hanya mengulang yang tersisa, itu akan memerlukan sedikit lebih banyak penyesuaian.

# Generate our custom methods
mice.impute.pmm_x1 <- 
  function (y, ry, x, donors = 5, type = 1, ridge = 1e-05, version = "", 
            ...) 
  {
    max_sum <- sum(max(x[,"x2"], na.rm=TRUE),
                   max(x[,"x3"], na.rm=TRUE))
    repeat{
      vals <- mice.impute.pmm(y, ry, x, donors = 5, type = 1, ridge = 1e-05,
                              version = "", ...)
      if (all(vals < max_sum)){
        break
      }
    }
    return(vals)
  }

mice.impute.pmm_x2 <- 
  function (y, ry, x, donors = 5, type = 1, ridge = 1e-05, version = "", 
            ...) 
  {
    repeat{
      vals <- mice.impute.pmm(y, ry, x, donors = 5, type = 1, ridge = 1e-05,
                              version = "", ...)
      if (all(vals == 0 | vals >= 14)){
        break
      }
    }
    return(vals)
  }

mice.impute.pmm_x3 <- 
  function (y, ry, x, donors = 5, type = 1, ridge = 1e-05, version = "", 
            ...) 
  {
    repeat{
      vals <- mice.impute.pmm(y, ry, x, donors = 5, type = 1, ridge = 1e-05,
                              version = "", ...)
      if (all(vals == 0 | vals >= 16)){
        break
      }
    }
    return(vals)
  }

Setelah kita selesai mendefinisikan metode, kita dengan mudah mengubah metode sebelumnya. Jika Anda hanya ingin mengubah satu variabel maka Anda cukup menggunakan imp_base$method["x2"] <- "pmm_x2"tetapi untuk contoh ini kami akan mengubah semua (penamaan tidak perlu):

imp_base$method <- c(x1 = "pmm_x1", x2 = "pmm_x2", x3 = "pmm_x3")

# The predictor matrix is not really necessary for this example
# but I use it just to illustrate in case you would like to 
# modify it
imp_ds <- 
  mice(dat, 
       method = imp_base$method, 
       predictorMatrix = imp_base$predictorMatrix)

Sekarang mari kita lihat dataset imputed ketiga:

> complete(imp_ds, action = 3)
    x1 x2 x3
1   21  0  0
2   50 19  0
3   31 18  0
4   15  0  0
5   36 19  0
6   82  0 54
7   14  0  0
8   14  0  0
9   19  0  0
10  18  0  0
11  16  0  0
12  36  0  0
13 583  0  0
14  50 22  0
15  52 19  0
16  14  0  0
17  50 22  0
18  52  0  0
19  26  0  0
20  24  0  0

Ok, itu berhasil. Saya suka solusi ini karena Anda dapat membonceng di atas fungsi utama dan hanya menambahkan batasan yang menurut Anda bermakna.

Memperbarui

Untuk menegakkan pengekangan ketat @ t0x1n yang disebutkan dalam komentar, kami mungkin ingin menambahkan kemampuan berikut ke fungsi pembungkus:

1. Simpan nilai yang valid selama loop sehingga data dari yang sebelumnya, sebagian yang berhasil berjalan tidak dibuang
2. Mekanisme melarikan diri untuk menghindari loop tak terbatas
3. Mengembang kumpulan donor setelah mencoba x kali tanpa menemukan kecocokan yang cocok (ini terutama berlaku untuk pmm)

Ini menghasilkan fungsi pembungkus yang sedikit lebih rumit:

mice.impute.pmm_x1_adv <-   function (y, ry, 
                                      x, donors = 5, 
                                      type = 1, ridge = 1e-05, 
                                      version = "", ...) {
  # The mice:::remove.lindep may remove the parts required for
  # the test - in those cases we should escape the test
  if (!all(c("x2", "x3") %in% colnames(x))){
    warning("Could not enforce pmm_x1 due to missing column(s):",
            c("x2", "x3")[!c("x2", "x3") %in% colnames(x)])
    return(mice.impute.pmm(y, ry, x, donors = 5, type = 1, ridge = 1e-05,
                           version = "", ...))
  }

  # Select those missing
  max_vals <- rowSums(x[!ry, c("x2", "x3")])

  # We will keep saving the valid values in the valid_vals
  valid_vals <- rep(NA, length.out = sum(!ry))
  # We need a counter in order to avoid an eternal loop
  # and for inflating the donor pool if no match is found
  cntr <- 0
  repeat{
    # We should be prepared to increase the donor pool, otherwise
    # the criteria may become imposs
    donor_inflation <- floor(cntr/10)
    vals <- mice.impute.pmm(y, ry, x, 
                            donors = min(5 + donor_inflation, sum(ry)), 
                            type = 1, ridge = 1e-05,
                            version = "", ...)

    # Our criteria check
    correct <- vals < max_vals
    if (all(!is.na(valid_vals) |
              correct)){
      valid_vals[correct] <-
        vals[correct]
      break
    }else if (any(is.na(valid_vals) &
                    correct)){
      # Save the new valid values
      valid_vals[correct] <-
        vals[correct]
    }

    # An emergency exit to avoid endless loop
    cntr <- cntr + 1
    if (cntr > 200){
      warning("Could not completely enforce constraints for ",
              sum(is.na(valid_vals)),
              " out of ",
              length(valid_vals),
              " missing elements")
      if (all(is.na(valid_vals))){
        valid_vals <- vals
      }else{
        valid_vals[is.na(valid_vals)] <- 
          vals[is.na(valid_vals)]
      }
      break
    }
  }
  return(valid_vals)
}

Perhatikan bahwa ini tidak bekerja dengan baik, kemungkinan besar karena set data yang disarankan gagal kendala untuk semua kasus tanpa hilang. Saya perlu menambah panjang loop menjadi 400-500 bahkan sebelum mulai berlaku. Saya berasumsi bahwa ini tidak disengaja, imputasi Anda harus meniru bagaimana data aktual dihasilkan.

Optimasi

Argumen tersebut ryberisi nilai-nilai yang tidak hilang dan kami mungkin dapat mempercepat loop dengan menghapus elemen yang telah kami temukan imputasi yang memenuhi syarat, tetapi karena saya tidak terbiasa dengan fungsi-fungsi internal saya telah menahan diri dari ini.

Saya pikir hal yang paling penting ketika Anda memiliki kendala kuat yang membutuhkan waktu untuk mengisi penuh adalah untuk memparalelasikan imputasi Anda ( lihat jawaban saya di CrossValidated ). Sebagian besar memiliki komputer saat ini dengan 4-8 core dan R hanya menggunakan salah satunya secara default. Waktu dapat (hampir) diiris menjadi dua dengan menggandakan jumlah inti.

Parameter yang hilang saat imputasi

Mengenai masalah x2hilang pada saat imputasi - tikus sebenarnya tidak pernah memasukkan nilai yang hilang ke dalam x- data.frame. The tikus metode termasuk mengisi beberapa nilai acak di awal. Bagian rantai dari imputasi membatasi dampak dari nilai awal ini. Jika Anda melihat mice-fungsi Anda dapat menemukan ini sebelum panggilan imputasi ( mice:::sampler-fungsi):

...
if (method[j] != "") {
  for (i in 1:m) {
    if (nmis[j] < nrow(data)) {
      if (is.null(data.init)) {
        imp[[j]][, i] <- mice.impute.sample(y, 
                                            ry, ...)
      }
      else {
        imp[[j]][, i] <- data.init[!ry, j]
      }
    }
    else imp[[j]][, i] <- rnorm(nrow(data))
  }
}
...

The data.initdapat dipasok ke micefungsi dan mice.imput.sample adalah prosedur dasar sampling.

Mengunjungi urutan

Jika urutan kunjungan penting, Anda dapat menentukan urutan di mana mice-fungsi menjalankan imputasi. Default adalah dari 1:ncol(data)tetapi Anda dapat mengatur visitSequenceuntuk menjadi apa pun yang Anda suka.

Hal terdekat yang bisa saya temukan adalah penyertaan informasi Amelia sebelumnya. Lihat bab 4.7 dalam sketsa , khususnya 4.7.2:

Priors tingkat pengamatan

Para peneliti sering memiliki informasi sebelumnya tambahan tentang nilai data yang hilang berdasarkan penelitian sebelumnya, konsensus akademik, atau pengalaman pribadi. Amelia dapat memasukkan informasi ini untuk menghasilkan imputasi yang jauh lebih baik. Algoritma Amelia memungkinkan pengguna untuk memasukkan prior Bayesian informatif tentang sel-sel data yang hilang individu dan bukan parameter model yang lebih umum, banyak di antaranya memiliki sedikit makna langsung.

Penggabungan prior mengikuti analisis Bayesian dasar di mana imputasi ternyata menjadi rata-rata tertimbang dari imputasi berbasis model dan rata-rata sebelumnya, di mana bobot adalah fungsi dari kekuatan relatif dari data dan sebelumnya: ketika model memprediksi dengan sangat baik , imputasi akan menurunkan bobot sebelumnya, dan sebaliknya (Honaker dan King, 2010).

Prior tentang pengamatan individu harus menggambarkan keyakinan analis tentang distribusi sel data yang hilang. Ini dapat berupa bentuk mean dan deviasi standar atau interval condensi. Sebagai contoh, kita mungkin tahu bahwa tarif tahun 1986 di Thailand sekitar 40%, tetapi kami memiliki beberapa ketidakpastian mengenai nilai pastinya. Keyakinan kami sebelumnya tentang distribusi sel data yang hilang, kemudian, berpusat pada 40 dengan deviasi standar yang mencerminkan jumlah ketidakpastian yang kami miliki tentang keyakinan kami sebelumnya.

Untuk memasukkan prior, Anda harus membuat matriks prior dengan empat atau lima kolom. Setiap baris matriks mewakili prior pada salah satu pengamatan atau satu variabel. Pada setiap baris, entri pada kolom pertama adalah baris pengamatan dan entri adalah kolom kedua adalah kolom observasi. Dalam empat matriks priors matrix kolom ketiga dan keempat adalah mean dan standar deviasi dari distribusi sebelumnya dari nilai yang hilang.

Jadi, meskipun Anda tidak akan dapat secara umum mengatakan sesuatu seperti x1<x2+x3, Anda bisa mengulang set data Anda dan menambahkan tingkat observasi sebelum untuk setiap kasus yang relevan. Batas konstan juga dapat diterapkan (seperti pengaturan x1, x2, dan x3 menjadi non-negatif). Sebagai contoh:

priors = matrix(NA, nrow=0, ncol=5);
for (i in seq(1, length(data))) 
{
    x1 = data$x1[i];
    x2 = data$x2[i];
    x3 = data$x3[i];

    if (is.na(x1) && !is.na(x2) && !is.na(x3))
    {
        priors = rbind(priors, c(i, 1, 0, x2+x3, 0.999999))
    }
}

amelia(data, m=1, bound = rbind(c(1, 0, Inf), c(2, 0, Inf), c(3, 0, Inf)), pr = priors);

Kendala mungkin lebih mudah diimplementasikan dalam prediksi pencocokan rata-rata imputasi ganda Ini mengasumsikan bahwa ada sejumlah besar pengamatan dengan variabel kendala yang tidak hilang yang memenuhi kendala. Saya sedang berpikir tentang mengimplementasikan ini dalam fungsi Hmiscpaket R. aregImputeAnda mungkin ingin memeriksa kembali dalam sebulan atau lebih. Penting untuk menentukan jarak maksimum dari target yang bisa diamati oleh donor, karena kendala akan mendorong donor lebih jauh dari donor ideal yang tidak dibatasi.

Saya percaya bahwa paket Amelia(Amelia II) saat ini memiliki dukungan paling komprehensif untuk menentukan batasan rentang nilai data. Namun, masalahnya adalah Ameliaasumsi bahwa data normal multivarian.

Jika dalam kasus Anda asumsi normalitas multivarian tidak berlaku, Anda mungkin ingin memeriksa micepaket, yang mengimplementasikan beberapa imputasi (MI) melalui persamaan dirantai . Paket ini tidak memiliki asumsi normalitas multivarian . Ini juga memiliki fungsi yang mungkin cukup untuk menentukan batasan , tapi saya tidak yakin sampai sejauh mana. Fungsi ini disebut squeeze(). Anda dapat membacanya di dokumentasi: http://cran.r-project.org/web/packages/mice/mice.pdf . Manfaat tambahan dari micefleksibilitasnya adalah dalam hal memungkinkan spesifikasi fungsi imputasi yang ditentukan pengguna dan pemilihan algoritma yang lebih luas. Berikut tutorial tentang melakukan MI, menggunakan mice:http://www.ats.ucla.edu/stat/r/faq/R_pmm_mi.htm .

Sejauh yang saya mengerti, Hmiscpaket Dr. Harrell , menggunakan pendekatan persamaan dirantai yang sama ( prediktif pencocokan rata-rata ), mungkin mendukung data yang tidak normal (dengan pengecualian normpmmmetode). Mungkin dia sudah menerapkan fungsionalitas spesifikasi kendala sesuai jawabannya di atas. Saya belum pernah menggunakan aregImpute(), jadi tidak bisa mengatakan lebih banyak tentang itu (saya sudah menggunakan Ameliadan mice, tapi saya jelas bukan ahli dalam statistik, hanya mencoba belajar sebanyak yang saya bisa).

Akhirnya, Anda mungkin menemukan ikhtisar berikut yang menarik, sedikit tanggal, tetapi masih bagus, tentang pendekatan, metode, dan perangkat lunak untuk beberapa imputasi data dengan nilai yang hilang: http://www.ncbi.nlm.nih.gov/pmc/articles / PMC1839993 . Saya yakin bahwa ada makalah tinjauan umum yang lebih baru tentang MI, tetapi hanya itu yang saya sadari saat ini. Saya harap ini agak membantu.

Jika saya memahami pertanyaan Anda dengan benar, menurut saya Anda sudah tahu nilai apa yang harus diambil oleh variabel yang hilang dengan tunduk pada beberapa kendala. Saya tidak terlalu fasih dalam SPSS tetapi di RI pikir Anda dapat menulis fungsi untuk melakukan itu (yang seharusnya tidak terlalu sulit tergantung pada pengalaman Anda, saya harus mengatakan). Saya tidak tahu ada paket yang berfungsi dengan kendala seperti itu.