---
title: "Test_factorization_methods_6"
output: html_document
date: "2024-11-06"
---

```{r setup, include=FALSE}
knitr::opts_chunk$set(echo = TRUE)
```


```{r}
library(edgeR)        # Version 3.38.1
library(fastTopics)   # Poisson NMF
library(NMF)          # Standard NMF
library(fcmclustering) # Fuzzy C-means clustering
library(mclust)       # Gaussian mixture models
library(WGCNA)        # Weighted gene co-expression network analysis
library(fcmclustering)
```






```{r}
# Read in counts matrix
TCGA_counts <- read.csv("../Data/TMM_TCGA_counts.csv", row.names=1)
CPTAC_counts <- read.csv("../Data/TMM_CPTAC_counts.csv", row.names=1)


# Log counts + 1 transform data for downstream use
log_TCGA_counts = log10(TCGA_counts + 1)
log_CPTAC_counts = log10(CPTAC_counts + 1)
```


```{r}
seed = 1 
set.seed(1)

number = 12

fit <- fit_poisson_nmf(as.matrix(TCGA_counts),k = number,numiter = 100,verbose = "progressbar",
                       control = list(extrapolate = TRUE))
fit <- poisson2multinom(fit)

# Write TCGA factorization with seed information included
write.table(fit$L, file=paste0("../Data/fastTopics/pNMF_H_TCGA_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
write.table(fit$F, file=paste0("../Data/fastTopics/pNMF_W_TCGA_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
```



```{r}
fit <- fit_poisson_nmf(as.matrix(log_TCGA_counts),k = number,numiter = 100,verbose = "progressbar",
                       control = list(extrapolate = TRUE))
fit <- poisson2multinom(fit)

# Write TCGA factorization with seed information included
write.table(fit$L, file=paste0("../Data/fastTopics/pNMF_H_TCGA_logscsub_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
write.table(fit$F, file=paste0("../Data/fastTopics/pNMF_W_TCGA_logscsub_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
```





```{r}
set.seed(1)
number <- 12
nmf_result <- nmf(TCGA_counts, rank = number, method = "brunet")

# Save the W matrix (genes × factors)
W <- basis(nmf_result)
write.table(W, file = "../Data/NMF/NMF_W_TCGA_genes.csv", 
            row.names = TRUE, col.names = TRUE, sep = ",")

# Save the H matrix (factors × samples)
H <- coef(nmf_result)
write.table(H, file = "../Data/NMF/NMF_H_TCGA_samples.csv", 
            row.names = TRUE, col.names = TRUE, sep = ",")




set.seed(1)
number <- 12
nmf_result <- nmf(log_TCGA_counts, rank = number, method = "brunet")

# Save the W matrix (genes × factors)
W <- basis(nmf_result)
write.table(W, file = "../Data/NMF/log_NMF_W_TCGA_genes.csv", 
            row.names = TRUE, col.names = TRUE, sep = ",")

# Save the H matrix (factors × samples)
H <- coef(nmf_result)
write.table(H, file = "../Data/NMF/log_NMF_H_TCGA_samples.csv", 
            row.names = TRUE, col.names = TRUE, sep = ",")
```






```{r}
set.seed(1)
number <- 12
seed <- 1

# FCM clustering
fcm_result <- fcm(TCGA_counts, centers=number, m=2)

# Get sample memberships (samples × clusters)
sample_memberships <- fcm_result$u
write.table(sample_memberships, 
            file=paste0("../Data/FCM/FCM_H_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")

# Get gene contributions using cluster centers
gene_weights <- fcm_result$centers
write.table(gene_weights, 
            file=paste0("../Data/FCM/FCM_W_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")





set.seed(1)
number <- 12
seed <- 1

# FCM clustering
fcm_result <- fcm(log_TCGA_counts, centers=number, m=2)

# Get sample memberships (samples × clusters)
sample_memberships <- fcm_result$u
write.table(sample_memberships, 
            file=paste0("../Data/FCM/log_FCM_H_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")

# Get gene contributions using cluster centers
gene_weights <- fcm_result$centers
write.table(gene_weights, 
            file=paste0("../Data/FCM/log_FCM_W_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")
```














```{r}
library(mclust)
set.seed(1)
number <- 12
seed <- 1  # Define seed variable

fit_samples <- Mclust(as.matrix(TCGA_counts), G=number)
probs_samples <- fit_samples$z
write.table(probs_samples, file=paste0("../Data/GMM/GMM_H_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")


cluster_centers <- fit_samples$parameters$mean  # This is (m x 12) - genes × clusters

gene_weights <- scale(cluster_centers, center = FALSE, scale = colSums(abs(cluster_centers)))

write.table(gene_weights, file=paste0("../Data/GMM/GMM_W_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")



library(mclust)
set.seed(1)
number <- 12
seed <- 1  # Define seed variable

fit_samples <- Mclust(as.matrix(log_TCGA_counts), G=number)
probs_samples <- fit_samples$z
write.table(probs_samples, file=paste0("../Data/GMM/log_GMM_H_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")


cluster_centers <- fit_samples$parameters$mean  # This is (m x 12) - genes × clusters

gene_weights <- scale(cluster_centers, center = FALSE, scale = colSums(abs(cluster_centers)))

write.table(gene_weights, file=paste0("../Data/GMM/log_GMM_W_TCGA_",toString(number),"_",toString(seed),".csv"), 
            row.names = TRUE, col.names = TRUE, sep = ",")
```


```{r}
set.seed(1)
number = 12

# Fit GMM
fit <- Mclust(as.matrix(log_TCGA_counts), G=number)
# Get cluster probabilities
probs <- fit$z

# Write TCGA factorization with seed information included
write.table(probs, file=paste0("../Data/GMM/GMM_W_TCGA_logscsub_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")

```






```{r}
library(WGCNA)

set.seed(1)
number = 12

# Create adjacency matrix
adjacency <- adjacency(as.matrix(TCGA_counts), power=number)
# Convert to topological overlap
TOM <- TOMsimilarity(adjacency)
# Identify modules
modules <- dynamicTreeCut(hierTOM)
```

```{r}
# Write TCGA factorization with seed information included
write.table(result@sampleFactors, file=paste0("../Data/WGCNA/WGCNA_H_TCGA_logscsub_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
write.table(result@featureLoadings, file=paste0("../Data/WGCNA/WGCNA_W_TCGA_logscsub_",toString(number),"_",toString(seed),".csv"), row.names = TRUE, col.names = TRUE, sep = ",")
```






