Hello! I am trying to perform a GLM FIR analysis and am having an issue with plotting data. I found the tutorial posted on the mne.tools page and based my analysis on that but am having issues with the data actually plotting. Everything runs well, I’m able to generate and see the actual data in a data frame, and all the variables I’m using are assigned the correct data/are not empty so I’m not sure why I’m having an issue. When I run the tutorial the data for that shows up on the plot so I’m not sure what’s going on. I’d appreciate any input, thanks so much!!
- MNE version: e.g. 1.6.1
- operating system: Windows 10
df["isR5N0Turn1"] = ["R5N0 Turn 1" in n for n in df["Condition"]]
df["isDelay"] = ["delay" in n for n in df["Condition"]]
df = df.query("isDelay in [True]")
df = df.query("isR5N0Turn1 in [True]")
# Make a new column that stores the condition name for tidier model below
# TidyCond is a way for us to keep things neat and organized
df.loc[:, "TidyCond"] = ""
df.loc[df["isR5N0Turn1"] == True, "TidyCond"] = "R5N0 Turn 1"
# Finally, extract the FIR delay in to its own column in data frame
df.loc[:, "delay"] = [n.split('_')[-1] for n in df.Condition]
# Only looking at our single condition
dm_cols_not_left = np.where(["R5N0 Turn 1" in c for c in dm.columns])[0]
dm = dm[[dm.columns[i] for i in dm_cols_not_left]]
lme = smf.mixedlm('theta ~ -1 + delay:TidyCond:Chroma', df,
groups=df["ID"]).fit()
fig, axes = plt.subplots(nrows=1, ncols=3, figsize=(20, 10))
# Extract design matrix columns that correspond to the condition of interest
dm_cond_idxs = np.where(["R5N0 Turn 1" in n for n in dm.columns])[0]
dm_cond = dm[[dm.columns[i] for i in dm_cond_idxs]]
# Extract the corresponding estimates from the lme dataframe for hbo
df_hbo = df_sum.query("TidyCond in ['R5N0 Turn 1']").query("Chroma in ['hbo']")
vals_hbo = [float(v) for v in df_hbo["Coef."]]
dm_cond_scaled_hbo = dm_cond * vals_hbo
# Extract the corresponding estimates from the lme dataframe for hbr
df_hbr = df_sum.query("TidyCond in ['R5N0 Turn 1']").query("Chroma in ['hbr']")
vals_hbr = [float(v) for v in df_hbr["Coef."]]
dm_cond_scaled_hbr = dm_cond * vals_hbr
# Extract the time scale for plotting.
# Set time zero to be the onset of the condition.
index_values = dm_cond_scaled_hbo.index - np.ceil(raw.annotations.onset[0])
index_values = np.asarray(index_values)
axes[0].plot(index_values, np.asarray(dm_cond))
axes[1].plot(index_values, np.asarray(dm_cond_scaled_hbo))
axes[2].plot(index_values, np.sum(dm_cond_scaled_hbo, axis=1), 'r')
axes[2].plot(index_values, np.sum(dm_cond_scaled_hbr, axis=1), 'b')
for ax in range(3):
axes[ax].set_xlim(-5, 30)
axes[ax].set_xlabel("Time (s)")
axes[0].set_ylim(-0.5, 1.3)
axes[1].set_ylim(-3, 8)
axes[2].set_ylim(-3, 8)
axes[0].set_title("FIR Model (Unscaled by GLM estimates)")
axes[1].set_title("FIR Components (Scaled by R5N0 Turn 1 GLM Estimates)")
axes[2].set_title("Evoked Response (R5N0 Turn 1)")
axes[0].set_ylabel("FIR Model")
axes[1].set_ylabel("Oyxhaemoglobin (ΔμMol)")
axes[2].set_ylabel("Haemoglobin (ΔμMol)")
axes[2].legend(["Oyxhaemoglobin", "Deoyxhaemoglobin"])